# Using Phylogenomic Data to Explore the Effects of Relaxed Clocks and Calibration Strategies on Divergence Time Estimation: Primates as a Test Case

Using Phylogenomic Data to Explore the Effects of Relaxed Clocks and Calibration Strategies on... Abstract Primates have long been a test case for the development of phylogenetic methods for divergence time estimation. Despite a large number of studies, however, the timing of origination of crown Primates relative to the Cretaceous–Paleogene (K–Pg) boundary and the timing of diversification of the main crown groups remain controversial. Here, we analysed a data set of 372 taxa (367 Primates and 5 outgroups, 3.4 million aligned base pairs) that includes nine primate genomes. We systematically explore the effect of different interpretations of fossil calibrations and molecular clock models on primate divergence time estimates. We find that even small differences in the construction of fossil calibrations can have a noticeable impact on estimated divergence times, especially for the oldest nodes in the tree. Notably, choice of molecular rate model (autocorrelated or independently distributed rates) has an especially strong effect on estimated times, with the independent rates model producing considerably more ancient age estimates for the deeper nodes in the phylogeny. We implement thermodynamic integration, combined with Gaussian quadrature, in the program MCMCTree, and use it to calculate Bayes factors for clock models. Bayesian model selection indicates that the autocorrelated rates model fits the primate data substantially better, and we conclude that time estimates under this model should be preferred. We show that for eight core nodes in the phylogeny, uncertainty in time estimates is close to the theoretical limit imposed by fossil uncertainties. Thus, these estimates are unlikely to be improved by collecting additional molecular sequence data. All analyses place the origin of Primates close to the K–Pg boundary, either in the Cretaceous or straddling the boundary into the Palaeogene. Bayes factors, Bayesian analysis, fossil, molecular dating, phylogenomic analysis, Primates, relaxed clock Divergence time estimation is fundamentally important to every field of evolutionary analysis. Reliable estimates of the timing of speciation events across the Tree of Life allow us to correlate these events with both biotic and abiotic phenomena on geological and more recent timescales, thus illuminating those that are most closely associated with periods of rapid diversification (Zhou et al. 2012; Andersen et al. 2015), evolutionary stasis (Alfaro et al. 2009; Liu et al. 2014), or high rates of extinction (Prum et al. 2015; Zhang et al. 2016). The prospect of employing genomic information for discovering the geological age of the Primates began virtually simultaneously with the publication of Zuckerkandl and Pauling’s (1965) molecular clock hypothesis. Sarich and Wilson (1967) employed a strict clock interpretation of immunological distance data to hypothesize that humans and other African apes (i.e., chimp and gorilla) shared a common ancestor as recently as 5 million years ago. This was revolutionary at the time given the implications for the necessarily rapid evolution of bipedal locomotion in the hominin lineage, and accordingly, drew considerable attention from the anthropological community (Read and Lestrel 1970; Uzzell and Pilbeam 1971; Lovejoy et al. 1972; Radinsky 1978; Corruccini et al. 1980). Despite this interest, it wasn’t until the 1980s that the field of divergence time estimation assumed a relatively modern flavor. It was then that investigators began to apply statistical models to DNA sequence data for the purposes of branch length and divergence time estimation (e.g., Hasegawa et al. 1985). Remarkably, these studies first emerged at a time when the sister–lineage relationship of humans to chimps was considered highly controversial (e.g., Goodman et al. 1983)—a relationship that is now considered unequivocal. Though it was Zuckerkandl and Pauling’s (1965) “molecular clock” hypothesis that made divergence time estimation feasible, it became clear very soon thereafter that there are myriad violations to a uniform clock. Thus, in subsequent decades, increasingly sophisticated models have been developed for relaxing the assumptions of a strict clock (Thorne et al. 1998; Kishino et al. 2001; Thorne and Kishino 2002; Drummond et al. 2006; Lepage et al. 2007; Rannala and Yang 2007; Heath et al. 2012). These models can be loosely divided into two categories: autocorrelated models, wherein rates of evolution in daughter species are statistically distributed around the parental rates (Sanderson 1997; Thorne et al. 1998; Aris-Brosou and Yang 2002), and uncorrelated models, wherein each lineage on the tree is free to assume a fully independent rate (Drummond et al. 2006; Rannala and Yang 2007; Paradis 2013). A parallel challenge to divergence time analysis can be observed in the development of calibration strategies (Marshall 1990; Yang and Rannala 2006; Benton and Donoghue 2007; Marshall 2008; Dornburg et al. 2011). Given that branch lengths on a phylogeny are the product of rate and time, investigators must assume one to infer the other. The most typical method for breaking this impasse is to employ fossil data as calibrations for one or more nodes in a given phylogeny so that the ages of other nodes can be inferred using DNA sequence data. This places an enormous burden on both the correct placement and age assignment of the fossils. If they are misplaced (i.e., assigned to clades where they do not belong) or if their geological ages are misinterpreted, divergence time estimates for the entire tree can be severely compromised (Martin 1993). The uncertainty imposed by paleontological ambiguity has not been as widely appreciated as have been the uncertainties introduced by the finite amount of DNA sequence data, which with the “genomics revolution,” is becoming steadily less problematic. We have reached a state of the art where branch lengths can be estimated with very high precision. The combination of genome-scale data, sophisticated molecular evolutionary models, and ever-increasing computational power has brought us to this point. Advances in DNA sequencing technology are such that virtually every major clade has at least a few species represented by whole-genome sequences, and this trend is rapidly accelerating. Bayesian methods have been developed such that parameter space can be explored during MCMC estimation, and though violations of the molecular clock will continue to present problems, methods for measuring and accommodating rate variation across phylogenies are explicit and generalizable. And finally, the computing power to integrate this information is increasing steadily. But because of the confounding effect of nonindependence of rate and time, the expectation of a conventional Bayesian analysis—that infinite data will eventually overcome prior information and provide posterior distributions with certainty—cannot be met (dos Reis and Yang 2013; Zhu et al. 2015). Thus, the field at present is focused on developing a better understanding of the effects of relaxed clock model choice, and on the impacts of calibration points, both with regard to abundance and placement in the phylogeny. Furthermore, in addressing these challenges, it is an open question as to whether simulation studies or tests of empirical data will be more informative for our understanding of best practices. With regard to clock model choice, an empirical study of three independent data sets showed that autocorrelated models outperform uncorrelated models, though the same study found a “high sensitivity” to the relaxation model employed (Lepage et al. 2007), while another empirical study found, however, that an independent rates model was superior (Linder et al. 2011). Simulation studies have only recently been employed, finding that even with complete taxon sampling, rate autocorrelation is challenging to detect (Ho et al. 2015). This has led to the conclusion that rigorous model selection should be conducted on a case-by-case basis, utilizing a wide range of real data sets, and thus comprising a promising avenue for future research (Ho et al. 2005; Ho and Duchene 2014; Ho et al. 2015). With regard to fossil calibration strategies, simulation studies (e.g., Duchene et al. 2014) have thus far agreed with previous empirical studies in finding that multiple calibrations are fundamentally important to accurate age estimation (Soltis et al. 2002; Yoder and Yang 2004; Linder et al. 2005; Rutschmann et al. 2007; Marshall 2008). Duchene et al. (2014) noted that calibrations close to the root of the phylogeny are most powerful. They also found that a significant source of error in divergence time estimation relates to misspecification of the clock model, especially when there are few informative calibrations. We cannot stress enough how sensitive posterior time estimates are to fossil information in constructing the time prior (Inoue et al. 2010). For example, different fossil calibration strategies have led to substantially different estimates of divergence times for mammals (dos Reis et al. 2014a) and early animals (dos Reis et al. 2015), regardless of how precise the branch length estimates are (Warnock et al. 2015). Thus, the field has reached the stage wherein there is general agreement that the choice of clock model and calibration strategy are fundamentally important to the accuracy of resulting age estimates, and thus, the way forward will clearly involve both empirical and simulation approaches to the problem. Here, we hope to contribute to this progress by conducting an exploration of model choice and calibration strategy in a classic empirical system: the Primates. Despite the fact that it is a relatively small and biologically uniform clade, primates have been inordinately and repeatedly the subject of divergence time analysis, with the first studies appearing at the very outset of molecular clock studies (Sarich and Wilson 1967), up to phylogenomic studies encompassing a large set of primate species (Perelman et al. 2011; Springer et al. 2012). This is largely due, undoubtedly, to the fact that we ourselves are members of this clade and can thus be forgiven for a persistent curiosity about our ancestral history. Age estimates for major primate divergence events have varied broadly among different studies (see Table 1), though one result has been relatively constant throughout: primate origins have been typically shown to predate the Cretaceous–Paleogene (K–Pg) mass-extinction event. Table 1. Overview of estimates of divergence times in Primates (in millions of years ago) for selected studies Study Data/Analysis Primates Haplorrhini Anthropoidea Platyrrhini Catarrhini Homininae Strepsirrhini Lorisiformes Lemuriformes Sarich and Wilson (1967) Immunological distance/strict clock 5 Hasegawa et al. (1985) 896 bp mtDNA/strict clock 3.7 Koop et al. (1989) B-globin DNA sequences/strict clock 55 40 25 Purvis (1995) Super Tree Analysis 57.2 39.9 14.7 8.1 41.8 22.1 39.6 Kumar and Hedges (1998) 658 nuclear genes (analyzed individually)/strict clock 47.6 6.7 Goodman et al. (1998) 60–80 kbp globin gene region/local clocks 63 58 40 7 50 23 45 Yang and Yoder (2003) 2404 bp mtDNA/local clocks 57.6 69.9 38.9 64.8 Poux and Douzery (2004) 1278 bp nDNA/local clocks 56.7, 58.4 45.4, 46.7 13.8, 14.2 39.6, 40.7 Eizirik et al. (2004) 8182 bp nDNA/relaxed clock 77.2 43.6 59.6 Steiper and Young (2006) 59.8 kbp of genomic data/ relaxed clock w/autocorrelated rates 77.5 42.9 20.8 30.5 8.6 57.1 Bininda-Emonds et al. (2007) 51,089 bp mtDNA and nDNA/ad hoc relaxed clock 87.7 Chatterjee et al. (2009) 6138 bp mtDNA and 2157 bp nDNA/relaxed clock 63.7 42.8 26.6 23.4 10.7 51.6 37.5 46.2 Perelman et al. (2011) 34,927 bp DNA/relaxed clock w/independent rates 87.2 81.3 43.5 24.8 31.6 8.3 68.7 40.3 58.6 Springer et al. (2012) 69 nDNA genes; 10 mtDNA genes/relaxed clock w/autocorrelated rates 67.8 61.2 40.6 23.3 25.1 8.0 54.2 34.7 50.0 dos Reis et al. (2012) 14,644 genes w/ 20.6 Mpb/relaxed clock w/autocorrelated rates 68.2 65.0 37.4 26.4 10.4 55.1 35.6 49.3 Finstermeier et al. (2013) complete mtDNA genomes/relaxed clock w/independent rates 66.2 63.0 45.3 22.0 32.0 8.4 56.9 34.5 47.1 Pozzi et al. (2014) complete mtDNA genomes/relaxed clock w/autocorrelated rates 74.1 70.0 46.7 20.9 32.1 10.6 66.3 40.3 59.6 This study (Strategy A) 3.4 million bp/relaxed clock w/ autocorrelated rates 74.4 70.6 45.0 25.3 32.6 10.5 62.7 37.9 57.2 This study (Strategy B) 3.4 million bp/relaxed clock w/ independent rates 84.8 78.6 45.9 25.3 28.8 7.6 64.0 38.2 55.3 Study Data/Analysis Primates Haplorrhini Anthropoidea Platyrrhini Catarrhini Homininae Strepsirrhini Lorisiformes Lemuriformes Sarich and Wilson (1967) Immunological distance/strict clock 5 Hasegawa et al. (1985) 896 bp mtDNA/strict clock 3.7 Koop et al. (1989) B-globin DNA sequences/strict clock 55 40 25 Purvis (1995) Super Tree Analysis 57.2 39.9 14.7 8.1 41.8 22.1 39.6 Kumar and Hedges (1998) 658 nuclear genes (analyzed individually)/strict clock 47.6 6.7 Goodman et al. (1998) 60–80 kbp globin gene region/local clocks 63 58 40 7 50 23 45 Yang and Yoder (2003) 2404 bp mtDNA/local clocks 57.6 69.9 38.9 64.8 Poux and Douzery (2004) 1278 bp nDNA/local clocks 56.7, 58.4 45.4, 46.7 13.8, 14.2 39.6, 40.7 Eizirik et al. (2004) 8182 bp nDNA/relaxed clock 77.2 43.6 59.6 Steiper and Young (2006) 59.8 kbp of genomic data/ relaxed clock w/autocorrelated rates 77.5 42.9 20.8 30.5 8.6 57.1 Bininda-Emonds et al. (2007) 51,089 bp mtDNA and nDNA/ad hoc relaxed clock 87.7 Chatterjee et al. (2009) 6138 bp mtDNA and 2157 bp nDNA/relaxed clock 63.7 42.8 26.6 23.4 10.7 51.6 37.5 46.2 Perelman et al. (2011) 34,927 bp DNA/relaxed clock w/independent rates 87.2 81.3 43.5 24.8 31.6 8.3 68.7 40.3 58.6 Springer et al. (2012) 69 nDNA genes; 10 mtDNA genes/relaxed clock w/autocorrelated rates 67.8 61.2 40.6 23.3 25.1 8.0 54.2 34.7 50.0 dos Reis et al. (2012) 14,644 genes w/ 20.6 Mpb/relaxed clock w/autocorrelated rates 68.2 65.0 37.4 26.4 10.4 55.1 35.6 49.3 Finstermeier et al. (2013) complete mtDNA genomes/relaxed clock w/independent rates 66.2 63.0 45.3 22.0 32.0 8.4 56.9 34.5 47.1 Pozzi et al. (2014) complete mtDNA genomes/relaxed clock w/autocorrelated rates 74.1 70.0 46.7 20.9 32.1 10.6 66.3 40.3 59.6 This study (Strategy A) 3.4 million bp/relaxed clock w/ autocorrelated rates 74.4 70.6 45.0 25.3 32.6 10.5 62.7 37.9 57.2 This study (Strategy B) 3.4 million bp/relaxed clock w/ independent rates 84.8 78.6 45.9 25.3 28.8 7.6 64.0 38.2 55.3 Notes: All estimates are mean estimates; see original works for confidence/credible intervals; all taxonomic designations signify crown nodes. Primates indicates the divergence of Strepsirrhini and Haplorrhini; Haplorrhini indicates the divergence of Tarsiidae from Anthropoidea; Platyrrhini; Catarrhini; Homininae indicates the divergence of Gorilla from Pan$$+$$Homo; Strepsirrhini indicates the divergence of Lorisiformes and Lemuriformes; Lorisiformes indicates the divergence of Galagidae and Lorisidae; Lemuriformes indicates the divergence of Daubentonia (the aye–aye) from other Malagasy lemurs. Table 1. Overview of estimates of divergence times in Primates (in millions of years ago) for selected studies Study Data/Analysis Primates Haplorrhini Anthropoidea Platyrrhini Catarrhini Homininae Strepsirrhini Lorisiformes Lemuriformes Sarich and Wilson (1967) Immunological distance/strict clock 5 Hasegawa et al. (1985) 896 bp mtDNA/strict clock 3.7 Koop et al. (1989) B-globin DNA sequences/strict clock 55 40 25 Purvis (1995) Super Tree Analysis 57.2 39.9 14.7 8.1 41.8 22.1 39.6 Kumar and Hedges (1998) 658 nuclear genes (analyzed individually)/strict clock 47.6 6.7 Goodman et al. (1998) 60–80 kbp globin gene region/local clocks 63 58 40 7 50 23 45 Yang and Yoder (2003) 2404 bp mtDNA/local clocks 57.6 69.9 38.9 64.8 Poux and Douzery (2004) 1278 bp nDNA/local clocks 56.7, 58.4 45.4, 46.7 13.8, 14.2 39.6, 40.7 Eizirik et al. (2004) 8182 bp nDNA/relaxed clock 77.2 43.6 59.6 Steiper and Young (2006) 59.8 kbp of genomic data/ relaxed clock w/autocorrelated rates 77.5 42.9 20.8 30.5 8.6 57.1 Bininda-Emonds et al. (2007) 51,089 bp mtDNA and nDNA/ad hoc relaxed clock 87.7 Chatterjee et al. (2009) 6138 bp mtDNA and 2157 bp nDNA/relaxed clock 63.7 42.8 26.6 23.4 10.7 51.6 37.5 46.2 Perelman et al. (2011) 34,927 bp DNA/relaxed clock w/independent rates 87.2 81.3 43.5 24.8 31.6 8.3 68.7 40.3 58.6 Springer et al. (2012) 69 nDNA genes; 10 mtDNA genes/relaxed clock w/autocorrelated rates 67.8 61.2 40.6 23.3 25.1 8.0 54.2 34.7 50.0 dos Reis et al. (2012) 14,644 genes w/ 20.6 Mpb/relaxed clock w/autocorrelated rates 68.2 65.0 37.4 26.4 10.4 55.1 35.6 49.3 Finstermeier et al. (2013) complete mtDNA genomes/relaxed clock w/independent rates 66.2 63.0 45.3 22.0 32.0 8.4 56.9 34.5 47.1 Pozzi et al. (2014) complete mtDNA genomes/relaxed clock w/autocorrelated rates 74.1 70.0 46.7 20.9 32.1 10.6 66.3 40.3 59.6 This study (Strategy A) 3.4 million bp/relaxed clock w/ autocorrelated rates 74.4 70.6 45.0 25.3 32.6 10.5 62.7 37.9 57.2 This study (Strategy B) 3.4 million bp/relaxed clock w/ independent rates 84.8 78.6 45.9 25.3 28.8 7.6 64.0 38.2 55.3 Study Data/Analysis Primates Haplorrhini Anthropoidea Platyrrhini Catarrhini Homininae Strepsirrhini Lorisiformes Lemuriformes Sarich and Wilson (1967) Immunological distance/strict clock 5 Hasegawa et al. (1985) 896 bp mtDNA/strict clock 3.7 Koop et al. (1989) B-globin DNA sequences/strict clock 55 40 25 Purvis (1995) Super Tree Analysis 57.2 39.9 14.7 8.1 41.8 22.1 39.6 Kumar and Hedges (1998) 658 nuclear genes (analyzed individually)/strict clock 47.6 6.7 Goodman et al. (1998) 60–80 kbp globin gene region/local clocks 63 58 40 7 50 23 45 Yang and Yoder (2003) 2404 bp mtDNA/local clocks 57.6 69.9 38.9 64.8 Poux and Douzery (2004) 1278 bp nDNA/local clocks 56.7, 58.4 45.4, 46.7 13.8, 14.2 39.6, 40.7 Eizirik et al. (2004) 8182 bp nDNA/relaxed clock 77.2 43.6 59.6 Steiper and Young (2006) 59.8 kbp of genomic data/ relaxed clock w/autocorrelated rates 77.5 42.9 20.8 30.5 8.6 57.1 Bininda-Emonds et al. (2007) 51,089 bp mtDNA and nDNA/ad hoc relaxed clock 87.7 Chatterjee et al. (2009) 6138 bp mtDNA and 2157 bp nDNA/relaxed clock 63.7 42.8 26.6 23.4 10.7 51.6 37.5 46.2 Perelman et al. (2011) 34,927 bp DNA/relaxed clock w/independent rates 87.2 81.3 43.5 24.8 31.6 8.3 68.7 40.3 58.6 Springer et al. (2012) 69 nDNA genes; 10 mtDNA genes/relaxed clock w/autocorrelated rates 67.8 61.2 40.6 23.3 25.1 8.0 54.2 34.7 50.0 dos Reis et al. (2012) 14,644 genes w/ 20.6 Mpb/relaxed clock w/autocorrelated rates 68.2 65.0 37.4 26.4 10.4 55.1 35.6 49.3 Finstermeier et al. (2013) complete mtDNA genomes/relaxed clock w/independent rates 66.2 63.0 45.3 22.0 32.0 8.4 56.9 34.5 47.1 Pozzi et al. (2014) complete mtDNA genomes/relaxed clock w/autocorrelated rates 74.1 70.0 46.7 20.9 32.1 10.6 66.3 40.3 59.6 This study (Strategy A) 3.4 million bp/relaxed clock w/ autocorrelated rates 74.4 70.6 45.0 25.3 32.6 10.5 62.7 37.9 57.2 This study (Strategy B) 3.4 million bp/relaxed clock w/ independent rates 84.8 78.6 45.9 25.3 28.8 7.6 64.0 38.2 55.3 Notes: All estimates are mean estimates; see original works for confidence/credible intervals; all taxonomic designations signify crown nodes. Primates indicates the divergence of Strepsirrhini and Haplorrhini; Haplorrhini indicates the divergence of Tarsiidae from Anthropoidea; Platyrrhini; Catarrhini; Homininae indicates the divergence of Gorilla from Pan$$+$$Homo; Strepsirrhini indicates the divergence of Lorisiformes and Lemuriformes; Lorisiformes indicates the divergence of Galagidae and Lorisidae; Lemuriformes indicates the divergence of Daubentonia (the aye–aye) from other Malagasy lemurs. Our study explores the effects of an autocorrelated versus an uncorrelated rate model on age estimates, and also explores the consequences of two different interpretations of both the age and the placement of key fossils with the living primate radiation. We apply these two strategies to a large phylogenomic data set for Primates (372 species and 3.4 million aligned base pairs). Until very recently, reliable calculations of branch lengths and age estimates within an analysis of this magnitude would have been beyond the capacity of computational methods. We have tackled many of these challenges by deploying the sequential Bayesian procedure used by dos Reis et al. (2012) wherein the posterior age estimates derived from a small taxonomic sample with genome-scale data are then deployed as priors for a subsequent analysis with many species and a much-reduced nucleotide sample. This procedure reduces the computational cost of a typical combined data analysis. It also helps to alleviate the concerns with the “missing data” problem as, in our approach, the sequence likelihood is only calculated for the species present in the alignment (Yang and Rannala 2006). The molecular timeline for primate evolution that emerges from this study can be interpreted with confidence. The data set is sufficiently large to provide highly precise branch length estimates, and the methods used are robust in accommodating violations of the molecular clock. The comparison of the two calibration strategies reveals their impact on the results by giving different age estimates within the tree, though the variation in inferred ages is not extreme. Which ages are considered most accurate will depend in large part on the degree of confidence in the fossils and their placement. As an unanticipated result of the study, the difference in age estimates for the deepest nodes of the phylogeny differ markedly when comparing the molecular rate models, with Bayesian model selection supporting the autocorrelated model. As with previous studies over the past several decades, the ancestral primate lineage is hypothesized to have survived the great K–Pg mass extinction event. Methods Bayesian estimates of divergence times of Primates were obtained using a supermatrix of molecular data with 372 species and 3.44 million base pairs (Mbp), combined with 17 fossil calibrations. The matrix is the result of merging the 372-species and 61 thousand base-pairs (kpb) data set of Springer et al. (2012) with a 10-species subset of the genome-scale alignment of dos Reis et al. (2012). Bayesian analyses were done using the program MCMCTree (Yang 2007). We assessed the robustness of time estimates by varying the clock model (strict clock, independent rates, and correlated rates) and by obtaining estimates under two fossil calibration strategies. Note that time estimates were obtained in two steps: in the first step, estimates were obtained for the small phylogeny of 10 species with a long alignment (3.38 Mbp). The marginal posterior of times was then used to construct the time prior in the second step for the 372-species phylogeny with a shorter alignment (61 kbp). This approach is almost the same as analysing the fully concatenated alignment in one step (3.38 Mbp $$+$$ 0.061 Mbp), but is computationally less demanding. All alignments, tree topology and fossil calibrations are available in Supplementary Material available on Dryad at http://dx.doi.org/10.5061/dryad.c020q. Sequence Alignment and Tree Topology Springer et al. (2012) alignment. We retrieved the sequence alignment of Springer et al. (2012), which is an extended version of the alignment of Perelman et al. (2011). The alignment has 372 species (367 primates and 5 outgroup species) and 79 gene segments (69 nuclear and 10 mitochondrial). The composite lagomorph sequence (an outgroup) was removed. We added a scandentian species (Tupaia belangeri), because its complete genome is available in the alignment of dos Reis et al. (2012), and because it has the 10 gene segments from the mitochondrial genome available (accession NC_002521). Many of the nuclear gene segments in the alignment of Springer et al. (2012) were mixtures of introns, exons, and UTRs, with out-of-frame indels in some exons. We manually curated the exons and separated the coding and noncoding segments of the alignment. These adjustments were necessary to facilitate an informed partition-based analysis of the data. Our modified version of Springer’s alignment was thus divided into six partitions: (i) 1st and 2nd codon positions for mitochondrial genes; (ii) 3rd positions for mitochondrial genes; (iii) mitochondrial RNA genes; (iv) 1st and 2nd codon positions for nuclear genes; (v) 3rd positions for nuclear genes; and (vi) noncoding segments of nuclear genes (UTRs and introns). The concatenated alignment has 372 species and is 61,132 base pairs long (Table 2). Our partitioning into codon positions and coding versus noncoding sequences follows established recommendations (Shapiro et al. 2006; Yang and Rannala 2006; Nascimento et al. 2017). Table 2. Sequence alignment summary Alignment Partition$$^{\mathrm{a}}$$ Sites Species Missing data$$^{\mathrm{b}}$$ Springer et al. 1. mit 1st$$+$$2nd 4816 330 61.4% 2. mit 3rd 2408 330 61.4% 3. mit RNA 2169 220 45.7% 4. nuclear 1st$$+$$2nd 16,309 239 53.8% 5. nuclear 3rd 8156 239 53.8% 6. nuclear noncoding 27,274 220 46.3% Partitions 1–6 61,132 372 51.1% (68.8%) dos Reis et al. 7. nuclear 1st$$+$$2nd 2,253,316 10 0.0% 8. nuclear 3rd 1,126,658 10 0.0% Partitions 7–8 3,379,974 10 0.0% (97.3%) Total 3,441,106 372 (96.8%) Alignment Partition$$^{\mathrm{a}}$$ Sites Species Missing data$$^{\mathrm{b}}$$ Springer et al. 1. mit 1st$$+$$2nd 4816 330 61.4% 2. mit 3rd 2408 330 61.4% 3. mit RNA 2169 220 45.7% 4. nuclear 1st$$+$$2nd 16,309 239 53.8% 5. nuclear 3rd 8156 239 53.8% 6. nuclear noncoding 27,274 220 46.3% Partitions 1–6 61,132 372 51.1% (68.8%) dos Reis et al. 7. nuclear 1st$$+$$2nd 2,253,316 10 0.0% 8. nuclear 3rd 1,126,658 10 0.0% Partitions 7–8 3,379,974 10 0.0% (97.3%) Total 3,441,106 372 (96.8%) $$^{\mathrm{a}}$$For topology estimation with RAxML, we used seven partitions: partitions 1 to 3, then 4 and 7 as one partition, 5 and 8 as one partition, and partition 6 divided into two: UTRs and Introns. $$^{\mathrm{b}}$$Numbers in brackets are the % missing data for the RAxML analysis. Note that MCMCTree only uses the species present in a partition to calculate the likelihood for the partition. In RAxML missing species in a partition are represented as sequences of only gaps in the partition, and thus the amount of missing data is larger. Table 2. Sequence alignment summary Alignment Partition$$^{\mathrm{a}}$$ Sites Species Missing data$$^{\mathrm{b}}$$ Springer et al. 1. mit 1st$$+$$2nd 4816 330 61.4% 2. mit 3rd 2408 330 61.4% 3. mit RNA 2169 220 45.7% 4. nuclear 1st$$+$$2nd 16,309 239 53.8% 5. nuclear 3rd 8156 239 53.8% 6. nuclear noncoding 27,274 220 46.3% Partitions 1–6 61,132 372 51.1% (68.8%) dos Reis et al. 7. nuclear 1st$$+$$2nd 2,253,316 10 0.0% 8. nuclear 3rd 1,126,658 10 0.0% Partitions 7–8 3,379,974 10 0.0% (97.3%) Total 3,441,106 372 (96.8%) Alignment Partition$$^{\mathrm{a}}$$ Sites Species Missing data$$^{\mathrm{b}}$$ Springer et al. 1. mit 1st$$+$$2nd 4816 330 61.4% 2. mit 3rd 2408 330 61.4% 3. mit RNA 2169 220 45.7% 4. nuclear 1st$$+$$2nd 16,309 239 53.8% 5. nuclear 3rd 8156 239 53.8% 6. nuclear noncoding 27,274 220 46.3% Partitions 1–6 61,132 372 51.1% (68.8%) dos Reis et al. 7. nuclear 1st$$+$$2nd 2,253,316 10 0.0% 8. nuclear 3rd 1,126,658 10 0.0% Partitions 7–8 3,379,974 10 0.0% (97.3%) Total 3,441,106 372 (96.8%) $$^{\mathrm{a}}$$For topology estimation with RAxML, we used seven partitions: partitions 1 to 3, then 4 and 7 as one partition, 5 and 8 as one partition, and partition 6 divided into two: UTRs and Introns. $$^{\mathrm{b}}$$Numbers in brackets are the % missing data for the RAxML analysis. Note that MCMCTree only uses the species present in a partition to calculate the likelihood for the partition. In RAxML missing species in a partition are represented as sequences of only gaps in the partition, and thus the amount of missing data is larger. dos Reis et al. (2012) alignment. We retrieved the genome-scale sequence alignment of dos Reis et al. (2012) of 36 mammal species, from which we extracted the sequences for 9 primates and 1 scandentian. The dos Reis et al. (2012) alignment was prepared using the highly curated mammalian genomes available in Ensembl. Though three additional primate genomes have become available in this database in the time since the original alignment was prepared, it is unlikely that their inclusion would change our results. The nine species represented in our study provide comprehensive phylogenetic representation of all major nodes in the primate tree and represent each of the higher-level clades (Fig. 1, inset). The original alignment has 14,632 nuclear, protein-coding genes, from which we removed 43 genes that were already present in the Springer alignment and 1 gene that was extremely long. All columns in the alignment with ambiguous nucleotides were removed, though care was taken not to disrupt the reading frame of the aligned coding sequences. The alignment was divided into two partitions: (i) 1st and 2nd codon positions; and (ii) 3rd codon positions. The final alignment has 10 species and is 3,441,106 base pairs long (missing data 0%, Table 2). Figure 1. View largeDownload slide The timetree of Primates. Nodes are drawn at their posterior mean ages in millions of years ago (Ma) estimated under the autocorrelated-rates (AR) clock model and calibration strategy A. Filled dots indicate nodes calibrated with the posterior times from the 10-species tree (inset figure), and empty dots indicate nodes with fossil constraints in the 372-species tree. Horizontal bars and numbers in parenthesis represent the 95% posterior CI for the node ages. Numbers associated with branches are ML bootstrap support values of major clades. Figure 1. View largeDownload slide The timetree of Primates. Nodes are drawn at their posterior mean ages in millions of years ago (Ma) estimated under the autocorrelated-rates (AR) clock model and calibration strategy A. Filled dots indicate nodes calibrated with the posterior times from the 10-species tree (inset figure), and empty dots indicate nodes with fossil constraints in the 372-species tree. Horizontal bars and numbers in parenthesis represent the 95% posterior CI for the node ages. Numbers associated with branches are ML bootstrap support values of major clades. Tree topology. The topology of the 372-species phylogeny was estimated by maximum likelihood (ML) using RAxML v 8.0.19 (Stamatakis 2014) under the GTR$$+$$G model (Yang 1994a,b), using seven partitions (Table 2) and 100 bootstrap replicates. Fossil Calibrations and Time Prior The two fossil calibration strategies used in this study are summarized in Table 3. They represent two different interpretations of the fossil record to construct calibrations for use in molecular clock dating analyses. Calibration strategy A is novel to this study, and calibration strategy B is based on the primate calibrations of dos Reis et al. (2012). Detailed justifications for the novel calibrations are provided in Appendix 1. Table 3. Fossil calibrations used in this study Calibration strategy$$^{\mathrm{a}}$$ Crown group Minimum (Ma) Maximum (Ma) MCMCTree Calibration$$^{\mathrm{b}}$$ Strategy A Human-Chimp$$^{\mathrm{c}}$$ 7.5 (†Sahelanthropus) 10 (unlikely before stem gorilla †Chororapithecus) B(0.075, 0.10, 0.01, 0.20) Human-Gorilla$$^{\mathrm{c}}$$ 10 (†Chororapithecus) 13.2 (unlikely before stem hominid †Sivapithecus) B(0.10, 0.132, 0.01, 0.20) Hominidae$$^{\mathrm{c}}$$ 11.2 (†Sivapithecus) 28 (unlikely before stem hominoid †Kamoyapithecus) B(0.112, 0.28, 0.01, 0.10) Catarrhini$$^{\mathrm{c}}$$ 25 (†Kamoyapithecus) 33.7 (absence of hominoids) B(0.25, 0.337, 0.01, 0.10) Anthropoidea$$^{\mathrm{c}}$$ 41 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 62.1 (K–Pg fossil modeling$$^{\mathrm{d}})$$ ST(0.4754, 0.0632, 0.98, 22.85) Strepsirrhini$$^{\mathrm{c}}$$ 37 (†Saharagalago) 58 (unlikely before †Altiatlasius) B(0.37, 0.58, 0.01, 0.10) Primates$$^{\mathrm{c}}$$ 57.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 88.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ S2N(0.474, 0.65, 0.0365, $$-$$3400, 0.6502, 0.1375, 11409) Euarchonta$$^{\mathrm{c}}$$ 65 (†Purgatorius) 130 (absence of placentals) G(36, 36.9) Strategy B Chimp-Human 5.7 (†Orrorin) 10 (absence of hominines) B(0.057, 0.10, 0.01, 0.05) Gorilla-Human 7.25 (†Chororapithecus) – L(0.0725, 0.1, 2) Hominidae 11.2 (†Sivapithecus) 33.7 (absence of pongines) B(0.112, 0.337, 0.05, 0.05) Catarrhini 23.5 (†Proconsul) 34 (absence of hominoids) B(0.235, 0.34, 0.01, 0.05) Anthropoidea 33.7 (†Catopithecus) — L(0.337, 0.1, 2) Haplorrrhini$$^{\mathrm{c}}$$ 45 (†Tarsius) — L(0.45, 0.1, 2) Strepsirrhini 33.7 (†Karanisia) 55.6 (absence of strepsirrhines) B(0.337, 0.556, 0.01, 0.05) Primates 55.6 (†Altiatlasius) — L(0.556, 0.1, 2) Euarchonta 61.5 (carpolestids and plesiadapids) 130 (absence of placentals) B(0.615, 1.30, 0.01, 0.05) Lorises$$^{\mathrm{c}}$$ 14 (Lorisidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.14, 0.37, 0.01, 0.10) Galagos$$^{\mathrm{c}}$$ 15 (Galagidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.15, 0.37, 0.01, 0.10) Lorisiformes$$^{\mathrm{c}}$$ 18 (†Mioeouticus) 38 (unlikely before †Karanisia/†Saharagalago) B(0.18, 0.38, 0.01, 0.10) Shared by Platyrrhini$$^{\mathrm{c}}$$ 15.7 (stem Pitheciinae) 33 (unlikely before †Catopithecus) B(0.157, 0.33, 0.01, 0.10) both Atelidae$$^{\mathrm{c}}$$ 12.8 (†Stirtonia) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) strategies Cebidae$$^{\mathrm{c}}$$ 12.8 (†Neosaimiri) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) Cercopithecinae$$^{\mathrm{c}}$$ 5 (†Parapapio) 23 (unlikely before †Prohylobates/Kamoyapithecus) B(0.05, 0.23, 0.01, 0.10) Colobinae$$^{\mathrm{c}}$$ 9.8 (†Microcolobus) 23 (unlikely before †Prohylobates/Kamoyapotheics) B(0.098, 0.23, 0.01, 0.10) Calibration strategy$$^{\mathrm{a}}$$ Crown group Minimum (Ma) Maximum (Ma) MCMCTree Calibration$$^{\mathrm{b}}$$ Strategy A Human-Chimp$$^{\mathrm{c}}$$ 7.5 (†Sahelanthropus) 10 (unlikely before stem gorilla †Chororapithecus) B(0.075, 0.10, 0.01, 0.20) Human-Gorilla$$^{\mathrm{c}}$$ 10 (†Chororapithecus) 13.2 (unlikely before stem hominid †Sivapithecus) B(0.10, 0.132, 0.01, 0.20) Hominidae$$^{\mathrm{c}}$$ 11.2 (†Sivapithecus) 28 (unlikely before stem hominoid †Kamoyapithecus) B(0.112, 0.28, 0.01, 0.10) Catarrhini$$^{\mathrm{c}}$$ 25 (†Kamoyapithecus) 33.7 (absence of hominoids) B(0.25, 0.337, 0.01, 0.10) Anthropoidea$$^{\mathrm{c}}$$ 41 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 62.1 (K–Pg fossil modeling$$^{\mathrm{d}})$$ ST(0.4754, 0.0632, 0.98, 22.85) Strepsirrhini$$^{\mathrm{c}}$$ 37 (†Saharagalago) 58 (unlikely before †Altiatlasius) B(0.37, 0.58, 0.01, 0.10) Primates$$^{\mathrm{c}}$$ 57.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 88.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ S2N(0.474, 0.65, 0.0365, $$-$$3400, 0.6502, 0.1375, 11409) Euarchonta$$^{\mathrm{c}}$$ 65 (†Purgatorius) 130 (absence of placentals) G(36, 36.9) Strategy B Chimp-Human 5.7 (†Orrorin) 10 (absence of hominines) B(0.057, 0.10, 0.01, 0.05) Gorilla-Human 7.25 (†Chororapithecus) – L(0.0725, 0.1, 2) Hominidae 11.2 (†Sivapithecus) 33.7 (absence of pongines) B(0.112, 0.337, 0.05, 0.05) Catarrhini 23.5 (†Proconsul) 34 (absence of hominoids) B(0.235, 0.34, 0.01, 0.05) Anthropoidea 33.7 (†Catopithecus) — L(0.337, 0.1, 2) Haplorrrhini$$^{\mathrm{c}}$$ 45 (†Tarsius) — L(0.45, 0.1, 2) Strepsirrhini 33.7 (†Karanisia) 55.6 (absence of strepsirrhines) B(0.337, 0.556, 0.01, 0.05) Primates 55.6 (†Altiatlasius) — L(0.556, 0.1, 2) Euarchonta 61.5 (carpolestids and plesiadapids) 130 (absence of placentals) B(0.615, 1.30, 0.01, 0.05) Lorises$$^{\mathrm{c}}$$ 14 (Lorisidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.14, 0.37, 0.01, 0.10) Galagos$$^{\mathrm{c}}$$ 15 (Galagidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.15, 0.37, 0.01, 0.10) Lorisiformes$$^{\mathrm{c}}$$ 18 (†Mioeouticus) 38 (unlikely before †Karanisia/†Saharagalago) B(0.18, 0.38, 0.01, 0.10) Shared by Platyrrhini$$^{\mathrm{c}}$$ 15.7 (stem Pitheciinae) 33 (unlikely before †Catopithecus) B(0.157, 0.33, 0.01, 0.10) both Atelidae$$^{\mathrm{c}}$$ 12.8 (†Stirtonia) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) strategies Cebidae$$^{\mathrm{c}}$$ 12.8 (†Neosaimiri) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) Cercopithecinae$$^{\mathrm{c}}$$ 5 (†Parapapio) 23 (unlikely before †Prohylobates/Kamoyapithecus) B(0.05, 0.23, 0.01, 0.10) Colobinae$$^{\mathrm{c}}$$ 9.8 (†Microcolobus) 23 (unlikely before †Prohylobates/Kamoyapotheics) B(0.098, 0.23, 0.01, 0.10) $$^{\mathrm{a}}$$Calibration strategies A and B are applied to the nodes in the phylogeny of 10 species. The shared calibrations are applied to the large tree of 372 species. Ages are in millions of years ago (Ma). $$^{\mathrm{b}}$$B($$t_{\mathrm{L}}, t_{\mathrm{U}}, p_{\mathrm{L}}, p_{\mathrm{U}})$$ means the node age is calibrated by a uniform distribution bounded between a minimum time $$t_{\mathrm{L}}$$, and a maximum time $$t_{\mathrm{U}}$$, with probabilities $$p_{\mathrm{L}}$$ and $$p_{\mathrm{U}}$$ that the age is outside the bounds. L($$t_{\mathrm{L}}$$, $$p$$, $$c$$, $$p_{\mathrm{L}})$$ means the node age is calibrated by a truncated Cauchy distribution with minimum age $$t_{\mathrm{L}}$$ and parameters $$p$$ and $$c$$, with the probability that the age is younger than the minimum bound to be $$p_{\mathrm{L}} \quad =$$ 5% (Inoue et al. 2010). ST($$a$$, $$b$$, $$c$$, $$d)$$ means the node age is calibrated by a skew-t density with parameters $$a$$, $$b$$, $$c,$$ and $$d$$ (Wilkinson et al. 2011). S2N($$a$$, $$b$$, $$c)$$ means the node age is calibrated by a mixture of two skew-normal distributions (Wilkinson et al. 2011). G($$a$$, $$b)$$ means the node age is calibrated by a gamma distribution with shape $$a$$ and rate $$b$$. $$^{\mathrm{c}}$$Detailed justifications for these fossil calibrations are given in the Appendix. For all other calibrations, justifications are in Benton et al. (2009) and dos Reis et al. (2012). $$^{\mathrm{d}}$$Calibration densities are the posterior distribution from a model of fossil preservation and discovery with species diversification that takes into account the effects of the K–Pg extinction in the model (Wilkinson et al. 2011). Table 3. Fossil calibrations used in this study Calibration strategy$$^{\mathrm{a}}$$ Crown group Minimum (Ma) Maximum (Ma) MCMCTree Calibration$$^{\mathrm{b}}$$ Strategy A Human-Chimp$$^{\mathrm{c}}$$ 7.5 (†Sahelanthropus) 10 (unlikely before stem gorilla †Chororapithecus) B(0.075, 0.10, 0.01, 0.20) Human-Gorilla$$^{\mathrm{c}}$$ 10 (†Chororapithecus) 13.2 (unlikely before stem hominid †Sivapithecus) B(0.10, 0.132, 0.01, 0.20) Hominidae$$^{\mathrm{c}}$$ 11.2 (†Sivapithecus) 28 (unlikely before stem hominoid †Kamoyapithecus) B(0.112, 0.28, 0.01, 0.10) Catarrhini$$^{\mathrm{c}}$$ 25 (†Kamoyapithecus) 33.7 (absence of hominoids) B(0.25, 0.337, 0.01, 0.10) Anthropoidea$$^{\mathrm{c}}$$ 41 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 62.1 (K–Pg fossil modeling$$^{\mathrm{d}})$$ ST(0.4754, 0.0632, 0.98, 22.85) Strepsirrhini$$^{\mathrm{c}}$$ 37 (†Saharagalago) 58 (unlikely before †Altiatlasius) B(0.37, 0.58, 0.01, 0.10) Primates$$^{\mathrm{c}}$$ 57.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 88.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ S2N(0.474, 0.65, 0.0365, $$-$$3400, 0.6502, 0.1375, 11409) Euarchonta$$^{\mathrm{c}}$$ 65 (†Purgatorius) 130 (absence of placentals) G(36, 36.9) Strategy B Chimp-Human 5.7 (†Orrorin) 10 (absence of hominines) B(0.057, 0.10, 0.01, 0.05) Gorilla-Human 7.25 (†Chororapithecus) – L(0.0725, 0.1, 2) Hominidae 11.2 (†Sivapithecus) 33.7 (absence of pongines) B(0.112, 0.337, 0.05, 0.05) Catarrhini 23.5 (†Proconsul) 34 (absence of hominoids) B(0.235, 0.34, 0.01, 0.05) Anthropoidea 33.7 (†Catopithecus) — L(0.337, 0.1, 2) Haplorrrhini$$^{\mathrm{c}}$$ 45 (†Tarsius) — L(0.45, 0.1, 2) Strepsirrhini 33.7 (†Karanisia) 55.6 (absence of strepsirrhines) B(0.337, 0.556, 0.01, 0.05) Primates 55.6 (†Altiatlasius) — L(0.556, 0.1, 2) Euarchonta 61.5 (carpolestids and plesiadapids) 130 (absence of placentals) B(0.615, 1.30, 0.01, 0.05) Lorises$$^{\mathrm{c}}$$ 14 (Lorisidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.14, 0.37, 0.01, 0.10) Galagos$$^{\mathrm{c}}$$ 15 (Galagidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.15, 0.37, 0.01, 0.10) Lorisiformes$$^{\mathrm{c}}$$ 18 (†Mioeouticus) 38 (unlikely before †Karanisia/†Saharagalago) B(0.18, 0.38, 0.01, 0.10) Shared by Platyrrhini$$^{\mathrm{c}}$$ 15.7 (stem Pitheciinae) 33 (unlikely before †Catopithecus) B(0.157, 0.33, 0.01, 0.10) both Atelidae$$^{\mathrm{c}}$$ 12.8 (†Stirtonia) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) strategies Cebidae$$^{\mathrm{c}}$$ 12.8 (†Neosaimiri) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) Cercopithecinae$$^{\mathrm{c}}$$ 5 (†Parapapio) 23 (unlikely before †Prohylobates/Kamoyapithecus) B(0.05, 0.23, 0.01, 0.10) Colobinae$$^{\mathrm{c}}$$ 9.8 (†Microcolobus) 23 (unlikely before †Prohylobates/Kamoyapotheics) B(0.098, 0.23, 0.01, 0.10) Calibration strategy$$^{\mathrm{a}}$$ Crown group Minimum (Ma) Maximum (Ma) MCMCTree Calibration$$^{\mathrm{b}}$$ Strategy A Human-Chimp$$^{\mathrm{c}}$$ 7.5 (†Sahelanthropus) 10 (unlikely before stem gorilla †Chororapithecus) B(0.075, 0.10, 0.01, 0.20) Human-Gorilla$$^{\mathrm{c}}$$ 10 (†Chororapithecus) 13.2 (unlikely before stem hominid †Sivapithecus) B(0.10, 0.132, 0.01, 0.20) Hominidae$$^{\mathrm{c}}$$ 11.2 (†Sivapithecus) 28 (unlikely before stem hominoid †Kamoyapithecus) B(0.112, 0.28, 0.01, 0.10) Catarrhini$$^{\mathrm{c}}$$ 25 (†Kamoyapithecus) 33.7 (absence of hominoids) B(0.25, 0.337, 0.01, 0.10) Anthropoidea$$^{\mathrm{c}}$$ 41 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 62.1 (K–Pg fossil modeling$$^{\mathrm{d}})$$ ST(0.4754, 0.0632, 0.98, 22.85) Strepsirrhini$$^{\mathrm{c}}$$ 37 (†Saharagalago) 58 (unlikely before †Altiatlasius) B(0.37, 0.58, 0.01, 0.10) Primates$$^{\mathrm{c}}$$ 57.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 88.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ S2N(0.474, 0.65, 0.0365, $$-$$3400, 0.6502, 0.1375, 11409) Euarchonta$$^{\mathrm{c}}$$ 65 (†Purgatorius) 130 (absence of placentals) G(36, 36.9) Strategy B Chimp-Human 5.7 (†Orrorin) 10 (absence of hominines) B(0.057, 0.10, 0.01, 0.05) Gorilla-Human 7.25 (†Chororapithecus) – L(0.0725, 0.1, 2) Hominidae 11.2 (†Sivapithecus) 33.7 (absence of pongines) B(0.112, 0.337, 0.05, 0.05) Catarrhini 23.5 (†Proconsul) 34 (absence of hominoids) B(0.235, 0.34, 0.01, 0.05) Anthropoidea 33.7 (†Catopithecus) — L(0.337, 0.1, 2) Haplorrrhini$$^{\mathrm{c}}$$ 45 (†Tarsius) — L(0.45, 0.1, 2) Strepsirrhini 33.7 (†Karanisia) 55.6 (absence of strepsirrhines) B(0.337, 0.556, 0.01, 0.05) Primates 55.6 (†Altiatlasius) — L(0.556, 0.1, 2) Euarchonta 61.5 (carpolestids and plesiadapids) 130 (absence of placentals) B(0.615, 1.30, 0.01, 0.05) Lorises$$^{\mathrm{c}}$$ 14 (Lorisidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.14, 0.37, 0.01, 0.10) Galagos$$^{\mathrm{c}}$$ 15 (Galagidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.15, 0.37, 0.01, 0.10) Lorisiformes$$^{\mathrm{c}}$$ 18 (†Mioeouticus) 38 (unlikely before †Karanisia/†Saharagalago) B(0.18, 0.38, 0.01, 0.10) Shared by Platyrrhini$$^{\mathrm{c}}$$ 15.7 (stem Pitheciinae) 33 (unlikely before †Catopithecus) B(0.157, 0.33, 0.01, 0.10) both Atelidae$$^{\mathrm{c}}$$ 12.8 (†Stirtonia) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) strategies Cebidae$$^{\mathrm{c}}$$ 12.8 (†Neosaimiri) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) Cercopithecinae$$^{\mathrm{c}}$$ 5 (†Parapapio) 23 (unlikely before †Prohylobates/Kamoyapithecus) B(0.05, 0.23, 0.01, 0.10) Colobinae$$^{\mathrm{c}}$$ 9.8 (†Microcolobus) 23 (unlikely before †Prohylobates/Kamoyapotheics) B(0.098, 0.23, 0.01, 0.10) $$^{\mathrm{a}}$$Calibration strategies A and B are applied to the nodes in the phylogeny of 10 species. The shared calibrations are applied to the large tree of 372 species. Ages are in millions of years ago (Ma). $$^{\mathrm{b}}$$B($$t_{\mathrm{L}}, t_{\mathrm{U}}, p_{\mathrm{L}}, p_{\mathrm{U}})$$ means the node age is calibrated by a uniform distribution bounded between a minimum time $$t_{\mathrm{L}}$$, and a maximum time $$t_{\mathrm{U}}$$, with probabilities $$p_{\mathrm{L}}$$ and $$p_{\mathrm{U}}$$ that the age is outside the bounds. L($$t_{\mathrm{L}}$$, $$p$$, $$c$$, $$p_{\mathrm{L}})$$ means the node age is calibrated by a truncated Cauchy distribution with minimum age $$t_{\mathrm{L}}$$ and parameters $$p$$ and $$c$$, with the probability that the age is younger than the minimum bound to be $$p_{\mathrm{L}} \quad =$$ 5% (Inoue et al. 2010). ST($$a$$, $$b$$, $$c$$, $$d)$$ means the node age is calibrated by a skew-t density with parameters $$a$$, $$b$$, $$c,$$ and $$d$$ (Wilkinson et al. 2011). S2N($$a$$, $$b$$, $$c)$$ means the node age is calibrated by a mixture of two skew-normal distributions (Wilkinson et al. 2011). G($$a$$, $$b)$$ means the node age is calibrated by a gamma distribution with shape $$a$$ and rate $$b$$. $$^{\mathrm{c}}$$Detailed justifications for these fossil calibrations are given in the Appendix. For all other calibrations, justifications are in Benton et al. (2009) and dos Reis et al. (2012). $$^{\mathrm{d}}$$Calibration densities are the posterior distribution from a model of fossil preservation and discovery with species diversification that takes into account the effects of the K–Pg extinction in the model (Wilkinson et al. 2011). Fossil calibration strategy A. We used the fossil-based prior densities constructed by Wilkinson et al. (2011) to calibrate the ages of crown Primates and crown Anthropoidea. The prior densities were constructed by modeling the processes of speciation, extinction, fossil preservation, and rates of fossil discovery in Primates. The effects of the K–Pg extinction were accounted for in the model. We calibrated six more node ages by using uniform distribution densities with soft bounds (Yang and Rannala 2006). We set the probability of violating a minimum bound to 1%. Because maximum bounds are based on weak evidence, we set the probability that a maximum bound is violated to 10% or 20%. The crown Haplorrhini node was left with no calibration as the branch separating that clade from crown Primates is very short and we wanted to avoid truncation with the fossil-modeling density on crown Primates. The prior on the age of crown Haplorrhini is instead set using the birth–death process with parameters $$\lambda =\mu =1$$ and $$\rho =0$$. These parameter values specify a uniform kernel density (Yang and Rannala 1997, equation 7). Fossil calibration strategy B. We used the same nine calibrations that dos Reis et al. (2012) used to calibrate the Primates and Scandentia clades. An additional calibration based on †Tarsius sp$$.$$ (Beard et al. 1994) was used for the Haplorrhini node. For nodes with a minimum bound only, modeled using a truncated Cauchy density, the spread parameter was set to $$c=2$$ (Inoue et al. 2010). For maximum bounds the probability that the bound was violated was set to 5%. There are other differences between strategies A and B (Table 1). For example, in A, we considered †Sahelanthropus, dated to 7.25 million years ago (Ma), to be the oldest member of the human-chimpanzee clade (Brunet et al. 2002, 2005) and used it to calibrate the clade accordingly, while in B, dos Reis et al. (2012) used †Orrorin (5.7 Ma) instead. In A, †Chororapithecus is given an age of 10 Ma (Geraads et al. 2002), while in B it is given the younger (perhaps more conservative) age of 7.25 Ma (Benton et al. 2009). We note that the ages of fossils and their relationships to extant groups are often controversial, and we cannot overemphasize the degree to which differences of opinion among palaeontologists are an important source of uncertainty in the construction of fossil calibrations, and accordingly, divergence time estimates throughout the phylogeny. Calibrating the 372-species phylogeny. Strategies A and B were used to obtain time estimates for the 10-species phylogeny using the 3.38 Mbp alignment. Then skew-t densities were fitted by ML to the marginal posterior ages of each of the 9 internal nodes in the 10-species phylogeny and used to calibrate the corresponding nodes in the 372-species tree. Eight additional fossil calibrations (Table 3) were used to calibrate additional nodes in the 372-species tree. For nodes without calibrations, the time prior was constructed using the birth–death process with parameters $$\lambda =\mu =1$$ and $$\rho =0$$. Bayesian time estimation then proceeded on the 372-species tree and 61 kbp alignment as usual. Rate Prior For the 10-species analysis, the rate prior was set as follows: the nuclear substitution rate at third codon positions in apes is roughly within 10$$^{\mathrm{-9}}$$ substitutions per site per year (s/s/y) (Burgess and Yang 2008). At first and second codon positions, it is about a fifth of the third position rate, or $$2\times 10^{\mathrm{-10}}$$ s/s/y. This gives roughly an overall rate of about $$5\times 10^{\mathrm{-10}}$$ s/s/y for the three positions combined. We thus used a diffuse gamma density G(2, 40) with mean 0.05 and 95% prior credibility interval (CI) 0.00606–0.139 (our time unit is 100 My, thus, this corresponds to $$6.06\times 10^{\mathrm{-12}}$$ to $$1.39\times 10^{\mathrm{-10}}$$ s/s/y). The analysis was conducted under both the autocorrelated rates (AR) and independent rates (IR) models. Parameter $$\sigma^{\mathrm{2}}$$ in the AR and IR models was assigned a gamma prior G(1, 10). Note that the average rate for loci, $$\mu_{i}$$, and $$\sigma_{i}^{2}$$ are assigned a gamma-Dirichlet prior (dos Reis et al. 2014b). For the 372-species phylogeny, the rate prior was assigned as follows: the mitochondrial substitution rate at third positions is about 20 times the rate at third positions in nuclear genes or $$2\times 10^{\mathrm{-8}}$$. Assuming 1st and 2nd codon positions evolve at about a fifth of the third position rate we get roughly $$4\times 10^{\mathrm{-9}}$$. The prior mean is then approximately $$2.5\times 10^{\mathrm{-9}}$$ s/s/y, which is the weighted average (by number of sites) of the substitution rates for the nuclear and mitochondrial partitions. We thus used a gamma density G(2, 8) with mean 0.25 and 95% CI 0.0302–0.696. For $$\sigma^{\mathrm{2}}$$ we used G(1, 10). MCMC and Bayesian Selection of Clock Model MCMC analyses were carried out with the program MCMCTree (Yang 2007), using the approximate likelihood method (dos Reis and Yang 2011). Convergence of the MCMC to the posterior distribution was assessed by running the analyses multiple times. MCMCtree runs were carried out without sequence data to calculate the joint prior of node ages. Results from all analyses were summarized as posterior means and 95% CIs. We have implemented marginal likelihood calculation by thermodynamic integration (path sampling) in the program MCMCTree. This allows us to calculate Bayes factors (BF) and posterior model probabilities to select for a clock model in the analysis. Details of our implementation are given in Appendix 2. Extensive discussions on marginal likelihood estimation by thermodynamic integration and stepping-stones (a related method) in the phylogenetic context are given in Lartillot and Philippe (2006), Lepage et al. (2007), Xie et al. (2011), and Baele et al. (2012a,b). A detailed simulation study is given in Ho et al. (2015). Thermodynamic integration is computationally intensive as we must sample from the power posterior $$f(\theta) f(D$$| $$\theta )^{\beta }$$, in a sampling path from the prior ($$\beta =0$$) to the posterior ($$\beta=1$$). Because the approximation to the likelihood is not good when samples are taken far away from the ML estimate (dos Reis and Yang 2011), as it happens when $$\beta$$ is small, the approximation cannot be used in the calculation of the power posterior. Thus, we use exact likelihood calculation on a smaller data set of nine primate species (inset of Fig. 1), for the six partitions of the Springer alignment (Table 2) to perform the Bayesian selection of clock model. We use 64 $$\quad$$$$\beta$$-points to construct the sampling path from the prior to the posterior and calculate the marginal likelihoods for the strict clock, and the AR and IR models. Effect of Genome-Scale Data In a conventional statistical inference problem, the variance of an estimate decreases in proportion to 1/$$n$$, with $$n$$ to be the sample size. Thus, as the sample size approaches infinity, the variance of an estimate approaches zero and the estimate converges to the true value. In divergence time estimation, which is an unconventional estimation problem, the nonidentifiability of times and rates means that the uncertainty in the posterior of times does not converge to zero as the amount of molecular data (the sample size) approaches infinity, but rather converge to a limiting value imposed by the uncertainties in the fossil calibrations (Yang and Rannala 2006; Rannala and Yang 2007). For infinitely long alignments, an infinite-sites plot (a plot of uncertainty in the time posterior, measured as the width of the CI, i.e., the difference between the 2.5% and 97.5% limits vs. the mean posterior of time) would converge onto a straight line. This line represents the amount of uncertainty in time estimates for every 1 million years (My) of divergence that is due solely to uncertainties in the fossil calibrations. We calculate the infinite-sites plots for time estimates on the 372-species phylogeny to study the effect of genome-scale data on the uncertainty of species-level time estimates. Results A Timeline of Primate Evolution Figure 1 illustrates time estimates under the AR model and calibration strategy A (Figs. 2–4 show detailed timetrees for the major clades). Under calibration strategy A and the AR model, we find that crown Primates originated 79.2–70.0 million years ago (Ma), before the K–Pg event at 66 Ma. However, the diversification of the main clades occurred much later. Crown Anthropoidea originated 48.3–41.8 Ma, with its two main crown groups, Catarrhini (Old World monkeys and apes) and Platyrrhini (New World monkeys) originating at 35.1–30.4 Ma and 27.5–23.6 Ma, respectively. Crown tarsiers originated 33.5–15.5 Ma. Crown Strepsirrhini date back to 66.8–58.8 Ma, with its two main crown groups, Lemuriformes and Lorisiformes, dating back to 61.6–52.7 Ma and 40.9–34.1 Ma, respectively. Figure 2. View largeDownload slide Strepsirrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 2. View largeDownload slide Strepsirrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 3. View largeDownload slide Catarrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 3. View largeDownload slide Catarrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 4. View largeDownload slide Tarsiidae and platyrrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 4. View largeDownload slide Tarsiidae and platyrrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Calibration strategy B under the AR model gives similar node age estimates for the younger nodes in the tree (i.e., the 95% CI of node age overlap, Fig. 5A). However, for the older nodes in the phylogeny (and in particular for Euarchonta, Primatomorpha, Primates, Haplorrhini, Lemuriformes, and Lemuriformes minus aye-aye), strategy A produced older estimates (Fig. 5A). Under strategy B a pre-K–Pg origin of crown Primates is also favored, although the posterior distribution of the age of crown Primates straddles the K–Pg boundary (71.4–63.9 Ma). The posterior probability for a pre-K–Pg origin of crown Primates is 80.0% under strategy B and 100% under strategy A. Posterior time estimates for all nodes under both strategies are given in the Supplementary Material spreadsheet available on Dryad. Figure 5. View largeDownload slide Effect of calibration strategy and relaxed-clock model. A) Posterior time estimates under fossil calibration strategy A versus time estimates under strategy B, for the AR clock model. B) Posterior time estimates under the AR versus IR clock models, for calibration strategy A. Figure 5. View largeDownload slide Effect of calibration strategy and relaxed-clock model. A) Posterior time estimates under fossil calibration strategy A versus time estimates under strategy B, for the AR clock model. B) Posterior time estimates under the AR versus IR clock models, for calibration strategy A. Note that the two calibration strategies are in many cases based on the same fossils (Table 3), and the intervals defined by the fossil bounds overlap extensively between the two strategies. However, the seemingly small differences between the two strategies lead to noticeable differences in the posterior time estimates (Fig. 5A). In general, minimum bound constraints are older in strategy A than in strategy B (Table 3), and thus this may be the cause of the older time estimates in A versus B. Time prior and effect of truncation and outgroups. User-specified calibration densities usually do not satisfy the constraint that descendant nodes must be younger than their ancestors, thus the dating methodology must “truncate” the calibration densities to satisfy the constraint to construct the time prior (Rannala 2016). The result is that user-specified calibration densities and marginal priors of node ages may look substantially different. Supplementary Figure S1 available on Dryad illustrates the effect of truncation on prior densities for strategies A and B. For example, in strategy B, the calibration densities on Euarchonta (the root of the phylogeny) and on Primates interact (the primate node has a Cauchy calibration density with a heavy tail), and consequently the prior density on the age of Euarchonta is pushed back (Supplementary Figure S1 available on Dryad). The result is that the marginal prior age of Euarchonta ranges from 136–178 Ma (Supplementary Figure S1 available on Dryad) instead of 130–61.5 Ma as in the calibration density (Table 3), while the upper age for the Primate prior is too old (127 Ma). In contrast, under strategy A, the calibration density on Primates has a much lighter tail, and thus the truncation effect with the Euarchonta node is minimal. The result is that the marginal time prior and the corresponding calibration densities for the Primates and Euarchonta nodes are very similar (Supplementary Figure S1 available on Dryad). Similarly, under strategy B, the priors for two other nodes (Anthropoidea and Human-Gorilla) that use the heavy-tailed Cauchy calibrations have upper 95% limits that also appear unreasonably old (86.3 Ma and 25.0 Ma, respectively). In general, calibration strategy A, which avoids using the long-tailed Cauchy calibrations, has calibration densities that are much closer to the resulting marginal priors, and thus strategy A results in a time prior which is much closer to the fossil information as interpreted by the palaeontologist. A set of calibrations for mitogenomic phylogenetic analysis. Mitogenomic markers are widely used to construct phylogenies of closely related primate species with examples seen in phylogeographic studies of diversification of primates in the Amazon (Nascimento et al. 2014) and in the timing of human diversification (Rieux et al. 2014). The posterior distributions obtained here for the 10-species genomic data are useful calibrations for mitogenomic studies. Note that these cannot be used if the molecular alignment contains nuclear data as the calibrations already contain the information from nuclear genomes. The list of skew-t calibrations is provided in Table 4, together with approximations based on the gamma distribution, which can be used in software that does not implement skew-t calibrations (such as BEAST or MrBayes). Table 4. Suggested skew-t and gamma calibrations for mitogenomic studies Crown group Skew-t$$^{\mathrm{a,b}}$$ Gamma$$^{\mathrm{a,c}}$$ Primates ST(0.878, 0.169, 2.41, 94.6) G(78.6, 77.6) Strepsirrhini ST(0.580, 0.0567, 2.12, 149) G(262, 327) Haplorrhini ST(0.705, 0.0628, 1.66, 960) G(254, 411) Anthropoidea ST(0.415, 0.0291, 0.949, 294) G(271, 363) Catarrhini ST(0.292, 0.0206, 0.995, 167) G(316, 733) Hominoidea ST(0.185, 0.0167, 2.44, 312) G(316, 1040) Human-Gorilla ST(0.0996, 0.0103, 19.6, 100) G(311, 1575) Human-Chimp ST(0.0788, 0.00687, 3.51, 6.15) G(292, 2715) Crown group Skew-t$$^{\mathrm{a,b}}$$ Gamma$$^{\mathrm{a,c}}$$ Primates ST(0.878, 0.169, 2.41, 94.6) G(78.6, 77.6) Strepsirrhini ST(0.580, 0.0567, 2.12, 149) G(262, 327) Haplorrhini ST(0.705, 0.0628, 1.66, 960) G(254, 411) Anthropoidea ST(0.415, 0.0291, 0.949, 294) G(271, 363) Catarrhini ST(0.292, 0.0206, 0.995, 167) G(316, 733) Hominoidea ST(0.185, 0.0167, 2.44, 312) G(316, 1040) Human-Gorilla ST(0.0996, 0.0103, 19.6, 100) G(311, 1575) Human-Chimp ST(0.0788, 0.00687, 3.51, 6.15) G(292, 2715) $$^{\mathrm{a}}$$Densities calculated under fossil calibration strategy A and AR model using a time unit of 100 My. $$^{\mathrm{b}}$$The parameters of the skew-t distribution are location, scale, shape, and df. $$^{\mathrm{c}}$$Note that here we use the shape ($$\alpha )$$ and rate ($$\beta )$$ parameterization. For the scale parameterization use $$s=1$$/$$\beta$$. The mean is $$\alpha$$/$$\beta$$ and variance is $$\alpha$$/$$\beta^{\mathrm{2}}$$. Table 4. Suggested skew-t and gamma calibrations for mitogenomic studies Crown group Skew-t$$^{\mathrm{a,b}}$$ Gamma$$^{\mathrm{a,c}}$$ Primates ST(0.878, 0.169, 2.41, 94.6) G(78.6, 77.6) Strepsirrhini ST(0.580, 0.0567, 2.12, 149) G(262, 327) Haplorrhini ST(0.705, 0.0628, 1.66, 960) G(254, 411) Anthropoidea ST(0.415, 0.0291, 0.949, 294) G(271, 363) Catarrhini ST(0.292, 0.0206, 0.995, 167) G(316, 733) Hominoidea ST(0.185, 0.0167, 2.44, 312) G(316, 1040) Human-Gorilla ST(0.0996, 0.0103, 19.6, 100) G(311, 1575) Human-Chimp ST(0.0788, 0.00687, 3.51, 6.15) G(292, 2715) Crown group Skew-t$$^{\mathrm{a,b}}$$ Gamma$$^{\mathrm{a,c}}$$ Primates ST(0.878, 0.169, 2.41, 94.6) G(78.6, 77.6) Strepsirrhini ST(0.580, 0.0567, 2.12, 149) G(262, 327) Haplorrhini ST(0.705, 0.0628, 1.66, 960) G(254, 411) Anthropoidea ST(0.415, 0.0291, 0.949, 294) G(271, 363) Catarrhini ST(0.292, 0.0206, 0.995, 167) G(316, 733) Hominoidea ST(0.185, 0.0167, 2.44, 312) G(316, 1040) Human-Gorilla ST(0.0996, 0.0103, 19.6, 100) G(311, 1575) Human-Chimp ST(0.0788, 0.00687, 3.51, 6.15) G(292, 2715) $$^{\mathrm{a}}$$Densities calculated under fossil calibration strategy A and AR model using a time unit of 100 My. $$^{\mathrm{b}}$$The parameters of the skew-t distribution are location, scale, shape, and df. $$^{\mathrm{c}}$$Note that here we use the shape ($$\alpha )$$ and rate ($$\beta )$$ parameterization. For the scale parameterization use $$s=1$$/$$\beta$$. The mean is $$\alpha$$/$$\beta$$ and variance is $$\alpha$$/$$\beta^{\mathrm{2}}$$. Effect of the Clock Model The clock model has a strong impact on posterior time estimates, particularly for the most ancient nodes in the phylogeny. Under the IR model, the ages of Euarchonta, Primatomorpha, Primates, Haplorrhini, and Strepsirrhini are substantially older than those estimated under the AR model (Fig. 5B). Posterior means and 95% CIs for locus (partition) rates obtained under both clock models are given in the Supplementary Material spreadsheet available on Dryad. Results of Bayesian model selection of clock model using thermodynamic integration are shown in Table 5. The AR model has the highest marginal likelihood in 5 of the 6 partitions analysed, with the posterior model probability > 90% in two partitions, and 79%, 66%, 53%, and 29% in the other four. When the six partitions are analysed in a multipartition data set, the posterior probability is virtually 100% in favor of the AR model. We note that ideally, the marginal likelihood calculations should have been carried out on the complete data set, but unfortunately, this is so computationally expensive that it cannot be done in a feasible amount of time. Thus, further work is necessary to confirm whether the preference for the AR model will remain in analysis of a more taxonomically dense primate phylogeny. Table 5. Bayesian selection of relaxed-clock model Data set Model Log marginal $$L^{\mathrm{a}}$$ BF$$^{\mathrm{b}}$$ $$P^{\mathrm{c}}$$ Mitochondrial 1st and 2nd c.p. SC –16,519.03 (0.010) $$1.3\times 10^{\mathrm{-18}}$$ $$1.2\times 10^{\mathrm{-18}}$$ IR –16,480.58 (0.021) 0.063 0.060 AR –16,477.82 (0.035) — 0.94 Mitochondrial 3rd c.p. SC –16,684.50 (0.014) — 0.61 IR –16,686.29 (0.043) 0.17 0.10 AR –16,685.26 (0.040) 0.47 0.29 Mitochondrial RNA SC –7906.85 (0.0087) 0.74 0.39 IR –7908.40 (0.015) 0.16 0.08 AR –7906.55 (0.023) — 0.53 Nuclear 1st and 2nd c.p. SC –32,179.80 (0.0092) 0.0047 0.0037 IR –32,175.77 (0.022) 0.27 0.21 AR –32,174.44 (0.032) — 0.79 Nuclear 3rd c.p. SC –24,535.33 (0.012) $$7.2\times 10^{\mathrm{-12}}$$ $$6.7\times 10^{\mathrm{-12}}$$ IR –24,512.45 (0.038) 0.062 0.058 AR –24,509.67 (0.030) — 0.94 Nuclear UTR and introns SC –64,739.20 (0.016) $$5.7\times 10^{\mathrm{-4}}$$ $$3.8\times 10^{\mathrm{-4}}$$ IR –64,732.41 (0.038) 0.51 0.34 AR –64,731.73 (0.046) — 0.66 All$$^{\mathrm{d}}$$ SC –162,684.8 (0.024) $$2.1 \times 10^{\mathrm{-103}}$$ $$2.1\times 10^{\mathrm{-103}}$$ IR –162,467.0 (0.086) $$8.4\times 10^{\mathrm{-9}}$$ $$8.4\times 10^{\mathrm{-9}}$$ AR –162,448.4 (0.15) — 1.00 Data set Model Log marginal $$L^{\mathrm{a}}$$ BF$$^{\mathrm{b}}$$ $$P^{\mathrm{c}}$$ Mitochondrial 1st and 2nd c.p. SC –16,519.03 (0.010) $$1.3\times 10^{\mathrm{-18}}$$ $$1.2\times 10^{\mathrm{-18}}$$ IR –16,480.58 (0.021) 0.063 0.060 AR –16,477.82 (0.035) — 0.94 Mitochondrial 3rd c.p. SC –16,684.50 (0.014) — 0.61 IR –16,686.29 (0.043) 0.17 0.10 AR –16,685.26 (0.040) 0.47 0.29 Mitochondrial RNA SC –7906.85 (0.0087) 0.74 0.39 IR –7908.40 (0.015) 0.16 0.08 AR –7906.55 (0.023) — 0.53 Nuclear 1st and 2nd c.p. SC –32,179.80 (0.0092) 0.0047 0.0037 IR –32,175.77 (0.022) 0.27 0.21 AR –32,174.44 (0.032) — 0.79 Nuclear 3rd c.p. SC –24,535.33 (0.012) $$7.2\times 10^{\mathrm{-12}}$$ $$6.7\times 10^{\mathrm{-12}}$$ IR –24,512.45 (0.038) 0.062 0.058 AR –24,509.67 (0.030) — 0.94 Nuclear UTR and introns SC –64,739.20 (0.016) $$5.7\times 10^{\mathrm{-4}}$$ $$3.8\times 10^{\mathrm{-4}}$$ IR –64,732.41 (0.038) 0.51 0.34 AR –64,731.73 (0.046) — 0.66 All$$^{\mathrm{d}}$$ SC –162,684.8 (0.024) $$2.1 \times 10^{\mathrm{-103}}$$ $$2.1\times 10^{\mathrm{-103}}$$ IR –162,467.0 (0.086) $$8.4\times 10^{\mathrm{-9}}$$ $$8.4\times 10^{\mathrm{-9}}$$ AR –162,448.4 (0.15) — 1.00 Notes: Marginal likelihoods are estimated by thermodynamic integration with 64 points. The substitution model is model is HKY$$+$$G. The age of the root is fixed to one (i.e., we use a “B(0.99, 1.01)” calibration on the root in MCMCTree). The rate priors are G(2, 1) and G(2, 20) for mitochondrial and nuclear data, respectively. The prior on $$\sigma^{\mathrm{2}}$$ is G(1, 1) in all cases. The model with the highest posterior probability in each data set is shown in bold type. SC $$=$$ strict clock; IR $$=$$ independent log-normal; AR $$=$$ autocorrelated rates.$$\sigma$$ $$^{\mathrm{a}}$$Values in brackets are the standard errors (see Appendix 2). $$^{\mathrm{b}}$$The values are the BF of the given model compared with the model with highest marginal likelihood (see Appendix 2). $$^{\mathrm{c}}$$Posterior model probabilities are calculated assuming a uniform prior on models. $$^{\mathrm{d}}$$The six data sets are analysed together as six partitions. Table 5. Bayesian selection of relaxed-clock model Data set Model Log marginal $$L^{\mathrm{a}}$$ BF$$^{\mathrm{b}}$$ $$P^{\mathrm{c}}$$ Mitochondrial 1st and 2nd c.p. SC –16,519.03 (0.010) $$1.3\times 10^{\mathrm{-18}}$$ $$1.2\times 10^{\mathrm{-18}}$$ IR –16,480.58 (0.021) 0.063 0.060 AR –16,477.82 (0.035) — 0.94 Mitochondrial 3rd c.p. SC –16,684.50 (0.014) — 0.61 IR –16,686.29 (0.043) 0.17 0.10 AR –16,685.26 (0.040) 0.47 0.29 Mitochondrial RNA SC –7906.85 (0.0087) 0.74 0.39 IR –7908.40 (0.015) 0.16 0.08 AR –7906.55 (0.023) — 0.53 Nuclear 1st and 2nd c.p. SC –32,179.80 (0.0092) 0.0047 0.0037 IR –32,175.77 (0.022) 0.27 0.21 AR –32,174.44 (0.032) — 0.79 Nuclear 3rd c.p. SC –24,535.33 (0.012) $$7.2\times 10^{\mathrm{-12}}$$ $$6.7\times 10^{\mathrm{-12}}$$ IR –24,512.45 (0.038) 0.062 0.058 AR –24,509.67 (0.030) — 0.94 Nuclear UTR and introns SC –64,739.20 (0.016) $$5.7\times 10^{\mathrm{-4}}$$ $$3.8\times 10^{\mathrm{-4}}$$ IR –64,732.41 (0.038) 0.51 0.34 AR –64,731.73 (0.046) — 0.66 All$$^{\mathrm{d}}$$ SC –162,684.8 (0.024) $$2.1 \times 10^{\mathrm{-103}}$$ $$2.1\times 10^{\mathrm{-103}}$$ IR –162,467.0 (0.086) $$8.4\times 10^{\mathrm{-9}}$$ $$8.4\times 10^{\mathrm{-9}}$$ AR –162,448.4 (0.15) — 1.00 Data set Model Log marginal $$L^{\mathrm{a}}$$ BF$$^{\mathrm{b}}$$ $$P^{\mathrm{c}}$$ Mitochondrial 1st and 2nd c.p. SC –16,519.03 (0.010) $$1.3\times 10^{\mathrm{-18}}$$ $$1.2\times 10^{\mathrm{-18}}$$ IR –16,480.58 (0.021) 0.063 0.060 AR –16,477.82 (0.035) — 0.94 Mitochondrial 3rd c.p. SC –16,684.50 (0.014) — 0.61 IR –16,686.29 (0.043) 0.17 0.10 AR –16,685.26 (0.040) 0.47 0.29 Mitochondrial RNA SC –7906.85 (0.0087) 0.74 0.39 IR –7908.40 (0.015) 0.16 0.08 AR –7906.55 (0.023) — 0.53 Nuclear 1st and 2nd c.p. SC –32,179.80 (0.0092) 0.0047 0.0037 IR –32,175.77 (0.022) 0.27 0.21 AR –32,174.44 (0.032) — 0.79 Nuclear 3rd c.p. SC –24,535.33 (0.012) $$7.2\times 10^{\mathrm{-12}}$$ $$6.7\times 10^{\mathrm{-12}}$$ IR –24,512.45 (0.038) 0.062 0.058 AR –24,509.67 (0.030) — 0.94 Nuclear UTR and introns SC –64,739.20 (0.016) $$5.7\times 10^{\mathrm{-4}}$$ $$3.8\times 10^{\mathrm{-4}}$$ IR –64,732.41 (0.038) 0.51 0.34 AR –64,731.73 (0.046) — 0.66 All$$^{\mathrm{d}}$$ SC –162,684.8 (0.024) $$2.1 \times 10^{\mathrm{-103}}$$ $$2.1\times 10^{\mathrm{-103}}$$ IR –162,467.0 (0.086) $$8.4\times 10^{\mathrm{-9}}$$ $$8.4\times 10^{\mathrm{-9}}$$ AR –162,448.4 (0.15) — 1.00 Notes: Marginal likelihoods are estimated by thermodynamic integration with 64 points. The substitution model is model is HKY$$+$$G. The age of the root is fixed to one (i.e., we use a “B(0.99, 1.01)” calibration on the root in MCMCTree). The rate priors are G(2, 1) and G(2, 20) for mitochondrial and nuclear data, respectively. The prior on $$\sigma^{\mathrm{2}}$$ is G(1, 1) in all cases. The model with the highest posterior probability in each data set is shown in bold type. SC $$=$$ strict clock; IR $$=$$ independent log-normal; AR $$=$$ autocorrelated rates.$$\sigma$$ $$^{\mathrm{a}}$$Values in brackets are the standard errors (see Appendix 2). $$^{\mathrm{b}}$$The values are the BF of the given model compared with the model with highest marginal likelihood (see Appendix 2). $$^{\mathrm{c}}$$Posterior model probabilities are calculated assuming a uniform prior on models. $$^{\mathrm{d}}$$The six data sets are analysed together as six partitions. Effect of Genome-Scale Data Figure 6 shows the infinite-sites plot for the primate data analysed here. For calibration strategy A, the eight primate nodes shared between the 10-species and 372-species trees (i.e., the nodes constrained by the large genome-scale alignment, Table 2) fall in an almost perfectly straight line ($$R \quad =$$ 0.992, Fig. 6A). This indicates that for these nodes, uncertainty in the time estimates is dominated by uncertainties in the fossil calibrations rather than by uncertainties in the molecular data. For strategy A, a regression line through the origin fitted to the eight data points is $$w \quad =$$ 0.128$$t$$, meaning that for every 1 My of divergence, 0.128 My are added to the CI width (Fig. 6A). On the other hand, when considering all 371 nodes in the tree, the relationship between CI width and mean times is far from linear, and the level of uncertainty is much higher. In this case, 0.277 My are added to the CI width for every 1 My of divergence. The trend is similar under calibration strategy B (Fig. 6B), albeit in this case there is in general more uncertainty in time estimates (i.e., the slope of the regression lines is larger). This appears due to strategy B being more conservative than strategy A, that is, some of the calibration densities used in B are substantially wider, encompassing larger time ranges (Supplementary Figure S1 available on Dryad). Figure 6. View largeDownload slide Infinite-sites plot. Posterior CI width is plotted against mean posterior divergence times for A) analysis under calibration strategy A and AR clock, and B) analysis under calibration strategy B and AR clock. In both cases, black dots indicate the eight primate nodes shared between the 10-species and 372-species trees, while the grey dots represent the rest of the nodes. Solid line: regression through the origin fitted to the black dots. Dashed line: regression through the origin fitted to all the dots. Figure 6. View largeDownload slide Infinite-sites plot. Posterior CI width is plotted against mean posterior divergence times for A) analysis under calibration strategy A and AR clock, and B) analysis under calibration strategy B and AR clock. In both cases, black dots indicate the eight primate nodes shared between the 10-species and 372-species trees, while the grey dots represent the rest of the nodes. Solid line: regression through the origin fitted to the black dots. Dashed line: regression through the origin fitted to all the dots. Discussion A Phylogenomic View of Primate Divergences A primary aim of this work was to study the effect of genome-scale data on divergence time estimates on a species-level phylogeny. Given the wide availability of whole genome data for a core-set of species, it is important to know whether the use of these data for a subsample of lineages will be enough to reduce time uncertainties in a species-level phylogeny to the theoretical limit. The results of Figure 6 clearly indicate this is not the case. Although for the core ancestral nodes in the primate phylogeny, the genome-scale alignments do constrain uncertainty in time estimates close to their theoretical limit (so it is highly unlikely that adding additional molecular data for these species will improve time estimates appreciably), for species without genome-scale data, there are still substantial uncertainties left for family-level and genus-level divergences in the tree. For some nodes, the CI-width is almost as large as the node age (e.g., for Tarsiidae, the node age is 23.8 Ma with CI-width 18 My, which is 76% of the node age, see Supplementary Material spreadsheet available on Dryad). Thus, much work is still needed in order to improve time estimates for the hundreds of more recent divergences in the tree. Furthermore, application of morphological-based models for dating (Ronquist et al. 2012) and the fossilized birth–death process (Heath et al. 2014) also offer exciting prospects and challenges in obtaining time estimates for the species-level divergences (O’Reilly et al. 2015; dos Reis et al. 2016). Improving these estimates will be important in studies of primate diversification rates and to correlate primate diversification events with major geological events in the history of the Earth (such as glaciations, continental drift, the closure the Panama isthmus, etc.). Sequential Bayesian Analysis Versus Secondary Calibrations In this work, we used the posterior of times obtained under a small data set as the prior of times in a second analysis under a large data set. This approach is justified as long as the data sets are independent under the likelihood model and as long as the data sets do not overlap (i.e., they share no genes). The use of the posterior in an analysis as the prior for the next is a well-known feature of Bayesian inference (Gelman et al. 2013). Consider data that can be split into two subsets, $$D = (D_{\mathrm{1}}$$, $$D_{\mathrm{2}})$$, which are independent under the likelihood model. The posterior distribution for parameter $$\theta$$ is ${\begin{array}{*{20}l} {f(\theta \vert D)} & \propto\, {f(\theta )f(D_{1} \vert \theta )f(D_{2} \vert \theta )}\\ & \propto\, {f(\theta \vert D_{1} )f(D_{2} \vert \theta )} \hfill \\ \end{array}},$ where $$f(\theta \vert D_{1} )\propto f(\theta )f(D_{1} \vert \theta )$$ is the posterior distribution of $$\theta$$ when only $$D_{\mathrm{1}}$$ are analysed. It is apparent that using $$f(\theta \vert D_{1} )$$ as the prior when analysing $$D_{\mathrm{2}}$$ leads to the posterior for the joint data $$D$$. In other words, performing the analysis in one step (joint analysis of $$D_{\mathrm{1}}$$ and $$D_{\mathrm{2}})$$ or in two steps (posterior under $$D_{\mathrm{1}}$$ as prior under $$D_{\mathrm{2}})$$ results in the same posterior distribution. The approach we used here to analyse the primate data is justified because the likelihood model assumes that the sequence partitions are nonoverlapping and independent. However, our approach is approximate. In multiparameter models, the posterior is a multidimensional distribution that may have a complex correlation structure. Here, we ignored the correlation structure of the nine times estimated using the genomic data and approximated the corresponding high-dimensional time posterior as the product of the marginal densities of the times, with a truncation applied to ensure that descendants are younger than ancestors. Note that joint analysis of all the partitions would have been preferable, but it is computationally prohibitive. This Bayesian sequential analysis is different from the use of secondary calibrations in some dating studies (Graur and Martin 2004), where the secondary calibrations were used as point calibrations (with the uncertainties on node age estimates ignored), and where in many cases the data analysed under the secondary calibration was the same as the data analysed to obtain the calibration in the first place (Graur and Martin 2004). Clock Model An interesting result from our study is the finding that the AR model fits the primate data better than the IR relaxed-clock model. In the context of previous studies, Lepage et al. (2007) found, using BF and no fossil calibrations, that two AR relaxed clocks (CIR and log-normal) fitted real data (eukaryotes, mammals, and vertebrates) better than IR models. More recently, Lartillot et al. (2016) introduced a mixed relaxed clock that has autocorrelated- and independent-rates components. In their analysis, the mixed clock appeared to provide a better description of rate evolution in the mammal phylogeny’ however, they did not assess clock model fit with BF. Linder et al. (2011) found, also by using BF, that IR models better fit an angiosperm phylogeny better than AR models. Additionally, they found that, when analysed without fossil calibrations, the AR model fit an ape phylogeny better than the IR model. However, when analysed with fossil calibrations, the IR model fit the ape data better. In the AR model the variance of the log-rate for branches is proportional to the time of divergence, so that the variance is expected to be close to zero for closely related species. In other words, the AR model allows for “local clocks” for closely related species, while allowing the rate to drift substantially across distantly related clades. This model is, from a biological point of view, quite appealing intuitively, and it also fits anecdotal evidence where the strict clock cannot be statistically rejected among very closely related species, for example, among the apes (dos Reis et al. 2016, box 2). In contrast, the IR model assumes that the variance of the branch rates is the same for different time scales. This would appear biologically unrealistic. Arguments have been put forward in favor and against both of the two types of relaxed-clock models examined by our study (Thorne et al. 1998; Drummond et al. 2006; Ho 2009), and clearly further research is still needed to understand the merits of each clock model and how they relate to evolutionary patterns in different genomic regions (Ho 2014). This will be a challenging task given how difficult it has been to distinguish between the two models in simulation studies (Heath et al. 2012; Ho et al. 2015; Lepage et al. 2007). Five Decades of Primate Molecular Timetrees Prior to the study reported here, the estimated age of the living primate clade has spanned a 30 My differential, ranging from as young as 55 Ma (Koop et al. 1989) to as old as 87 Ma (Perelman et al. 2011). The new millennium has been a particularly active time for primate divergence time analysis. Beginning in the early 2000’s, published studies have employed a diverse assortment of data sets applied to the problem (e.g., nuclear, mitochondrial, and their combination), as well as a range of statistical methods and calibration densities. Despite this array of data and methods, all of these studies—with only one notable outlier (Chatterjee et al. 2009)—have consistently indicated that the crown Primates clade originated prior to the K–Pg event (see also Steiper and Seiffert 2012). Given the continued dearth of fossil data to support this hypothesis, however, the result continues to be viewed with scepticism by the paleoanthropological community (Bloch et al. 2007; Silcox 2008; O’Leary et al. 2013; but see Martin et al. 2007). As described at length above, the current study gives added weight to the conclusion that primates are an ancient clade of placental mammals, arising just prior to or millions of years before the K–Pg. And even though lineage diversification within the major subclades is hypothesized not to have occurred until after the commencement of the Paleogene, the separation of tarsiers from other haplorrhines, and the divergence of haplorrhines and strepsirrhines, consistently appear to proceed or nearly coincide with the K–Pg. Given that this event was unequivocally one of the most disruptive and destructive geological episodes in Earth history, the temporal coincidence speaks both to the ecological flexibility and to the evolutionary opportunism of the earliest primates. Although now extinct in North America and in Europe, the primate fossil record shows that the clade was once nearly pan-global, even potentially including Antarctica. Thus, when viewed in the context of divergence date estimates, all of which fall within a temporal window when, as now, continental and island residences would have already been sundered by significant oceanic barriers (most notably, the separation of South America from Africa by the Atlantic Ocean), we must conclude that early primates would have been able dispersers. In fact, the ability to cross barriers, both terrestrial and aquatic, and to successfully colonize new land masses, are distinct hallmarks of the primate radiation (Gingerich 1981; Yoder and Nowak 2006; de Queiroz 2005, 2014; Seiffert 2012; Beard 2016; Bloch et al. 2016). Supplementary material Data available from the Dryad Digital Repository: http://dx.doi.org/10.5061/dryad.c020q. An R package, mcmc3r, that helps the user with marginal likelihood estimation with MCMCTree is available from https://github.com/dosreislab/mcmc3r. Funding This work was supported by grant BB/J009709/1 from the Biotechnology and Biosciences Research Council (UK). Part of this work was carried out while MdR was visiting the National Evolutionary Synthesis Center (NESCent, National Science Foundation #EF-0905606) in summer 2013. joint scholarship from University College London and the Government of Mexico’s CONACYT to J.B.M.; Grant from the Burroughs Wellcome Fund and by DEB-1354610 from the National Science Foundation to A.D.Y. Dedication This article is dedicated to the memory of our co-author Gregg F. Gunnell and to the lasting value of his work. Appendix 1 Justifications for dates assigned to 17 fossil calibrations in this study (Table 3) are given below: 8 calibrations for strategy A (SA), 1 calibration for strategy B (SB), and 8 calibrations shared among both strategies (SAB). The justifications for the remaining calibrations in Table 3 are given in dos Reis et al. (2012; see also Benton et al. 2009). In some cases, the dates used are not exactly those published in cited references. In these cases, the dates utilized reflect published as well as unpublished information or adjustments deemed necessary given the uncertainty of some dates. In any case, the discrepancies are always small and unimportant considering the breadth of the fossil calibrations. Note that specifying maximum bounds is a difficult task, in particular because absence of fossil evidence is not evidence that a clade did not exist in a point in time (Ho and Philips 2009). We use stem fossils as benchmark points onto which to construct diffuse maximum bounds (i.e., with a large probability of violation, $$p_{\mathrm{U}})$$ on some node ages. Hominini | Homo-Pan | 7.5 Ma–10 Ma | SA The minimum age for the divergence of hominins is placed at 7.5 Ma and is based on the appearance of †Sahelanthropus at 7.2 Ma (Brunet et al. 2002, 2005; Lebatard et al. 2008). There is some controversy as to the proper taxonomic position of Sahelanthropus (Wolpoff et al. 2006; MacLatchy et al. 2010), but we regard it as the oldest record of a plausible crown hominin. Sahelanthropus comes from the Anthracotheriid Unit of an unnamed formation in the Mega-Chad Basin in Chad (Brunet et al. 2005). The associated mammalian fauna is very similar to that found from the Nawata Formation, Lothagam, Kenya which may be as old as 7.4 Ma (MacDougall 2003). The divergence of the hominin lineage seems unlikely to have occurred before the appearance of the potential gorillin Chororapithecus at 10 Ma (Suwa et al. 2007; Harrison 2010a). Homininae | Gorilla-Homo | 10 Ma–13.2 Ma | SA The minimum age for the divergence of crown hominines is placed at 10 Ma based on the appearance of the potential gorillin Chororapithecus. Like Sahelanthropus, the taxonomic status of Chororapithecus is not without controversy (Suwa et al. 2007; Harrison 2010a). Harrison (2010a) regards Chororapithecus as best interpreted as a stem hominin or even a stem hominid—we feel that the features that do support a relationship with gorillas are well enough established to use the date of appearance of Chororapithecus as a minimum divergence date for hominines. Chororapithecus comes from the late Miocene Beticha section of the Chorora Formation in Ethiopia and is dated at 10–10.5 Ma (Geraads et al. 2002). We use a maximum divergence date of 13.2 (Raza et al. 1983) for Sivapithecus, but it is now evident that this date might be slightly too old. In light of this, the divergence of the hominine lineage is unlikely to have taken place before the earliest appearance of the probable crown pongine Sivapithecus at 12.5 Ma (Begun 2010; Begun et al. 2012). Hominidae | Pongo-Homo | 11.2 Ma–28 Ma | SA The minimum age for the divergence of crown hominids is placed at 11.2 Ma based on the earliest appearance of the crown pongine Sivapithecus (Kappelman et al. 1991; Begun 2010; Begun et al. 2012). Sivapithecus is known from Siwalik group rocks (Chinji, Nagri, and Dhok Pathan formations) in Indo-Pakistan that range in age from 14 Ma to 5.5 Ma (Badgley and Behrensmeyer 1995) with Sivapithecus restricted to a range of 12.5 Ma to 7.4 Ma (Flynn et al. 1995). The divergence of crown hominids is unlikely to have occurred before the first appearance of Kamoyapithecus at 25 Ma (Seiffert 2010). We use 28 Ma as a slightly more conservative maximum. Catarrhini | Homo-Macaca | 25 Ma–33.7 Ma | SA The presence of the crown hominoid Kamoyapithecus (Zalmout et al. 2010) indicates that the minimum divergence time for crown Catarrhini is 25 Ma (Seiffert 2010). Kamoyapithecus is only known from the Erageliet beds, Kalakol Basalts locality of Lothidok in Kenya (Madden 1980; Leakey et al. 1995; Rasmussen and Gutierrez 2009). A soft maximum of 33.7 Ma on the age of crown catarrhines is given due to the absence of hominoids before 33.7 Ma. Anthropoidea | Catarrhini–Platyrrhini | 41 Ma–62.1 Ma | SA Calibration density constructed from fossil modeling. The effects of the K–Pg extinction are included in the model. See Wilkinson et al. (2011) for details. Haplorrhini | crown Tarsius | 45 Ma | SB The presence of the crown tarsiid Xanthorhysis in Shanxi Province, China (Beard 1998) and apparently of the genus Tarsius in fissure fills at Shanghuang in Jiangsu Province, China (Beard et al. 1994), both dating to the late middle Eocene (40–45 Ma), circumscribe the minimum divergence time for crown Haplorrhini at 45 Ma. Strepsirrhini | Lorisiformes–Lemuriformes | 37 Ma–58 MA | SA The minimum age for the divergence of crown Strepsirrhini is placed at 37 Ma based on the first appearance of the crown lorisiform Saharagalago (Seiffert et al. 2003). Saharagalago is only known from Fayum Quarry BQ-2 in the Birket Qarun Formation, Egypt. The divergence of crown strepsirrhines is unlikely to have occurred before the first appearance of the basal primate Altiatlasius from Ouarzazate in Morocco which is considered to represent the late Paleocene (Thanetian) and dating to around 58 Ma (Gheerbrant et al. 1993, 1998; Seiffert 2010). Primates | Haplorrhini–Strepsirrhini | 57.6 Ma–88.6 Ma | SA Calibration density constructed from fossil modeling. Includes the effects of the K–Pg extinction in the model. See Wilkinson et al. (2011) for details. Euarchonta | Scandentia–Primates | 65 Ma–130 Ma | SA The minimum age for the divergence of crown Euarchonta is placed at 65 Ma based on the first appearance of the crown euarchontan Purgatorius (Bloch et al. 2007). Purgatorius is known from the early Paleocene (Puercan) Tullock and Bear Formations in Montana (Clemens 1974; Buckley 1997; Clemens 2004; Chester et al. 2015) and from the earliest Paleocene Ravenscrag Formation in Saskatchewan (Fox and Scott 2011). The divergence of Euarchonta is unlikely to have been before the appearance of placental mammals by at least 130 Ma (Luo 2007; but see Luo et al. 2011 for a potential 130 My old eutherian). Lorisidae | Nycticebus–Perodicticus | 14 Ma–37 Ma | SAB The minimum age for the divergence of crown Lorisidae is placed at 14 Ma based on an undescribed genus and species from Fort Ternan in Kenya cited by Walker (1978) and Harrison (2010b). The minimum age could possibly be as old as 19 Ma if Mioeuoticus (Leakey in Bishop 1962; Walker 1978) represents a crown lorisid (Harrison 2010b). Fossil lorisids are known from the early to middle Miocene in Africa (Phillips and Walker 2000, 2002; Harrison 2010b) and from the late Miocene of Pakistan (Jacobs 1981). The divergence of lorisoids is unlikely to have occurred before the first appearance of the potential stem lorisid Karanisia at BQ-2 in Egypt (Seiffert et al. 2003). Galagidae | Galago–Euoticus | 15 Ma–37 Ma | SAB The minimum age for the divergence of crown Galagidae is placed at 15 Ma based on an undescribed genus and species from Maboko Island in Kenya cited by McCrossin (1999) and Harrison (2010b). The minimum age could possibly be as old as 19 Ma if either Progalago or Komba represent a crown galagid (MacInnes 1943; Simpson 1967; Harrison 2010b). Fossil galagids are known from the early Miocene through early Pleistocene in Africa (Phillips and Walker 2002; Harrison 2010b). The divergence of galagids is unlikely to have occurred before the first appearance of the potential stem lorisid Karanisia at BQ-2 in Egypt (Seiffert et al. 2003). Lorisiformes | Nycticebus–Perodicticus | 14 Ma–37 Ma | SAB The minimum age for the divergence of crown Lorisiformes is placed at 18 Ma based on the appearance of the potential crown lorisoid Mioeuoticus in the early Miocene of East Africa (Harrison 2010b). The divergence of lorisiforms is unlikely to have occurred before the first appearances of Karanisia at BQ-2 in Egypt (Seiffert et al. 2003). We use 38 Ma as a conservative soft maximum. Platyrrhini | Pitheciidae–Callitrichidae | 15.7 Ma–33 Ma | SAB The minimum age for the divergence of crown Platyrrhini is based on the first occurrence of the crown pitheciine Proteropithecia dated at 15.7 Ma (Kay et al. 1998; Fleagle and Tejedor 2002). The minimum age could be as much 18 Ma if either (or both) Soriacebus or Carlocebus represent crown pitheciins (Fleagle et al. 1987, 1995; Fleagle 1990; Bown and Fleagle 1993; Rosenberger 2011). All of these taxa are known from the early and middle Miocene of Argentina (Fleagle and Tejedor 2002). The divergence of platyrrhines is unlikely to have occurred before the appearance of the crown catarrhine Catopithecus (33 Ma, Fayum, Egypt) although a recently published report has claimed a 36 Ma date for a stem platyrrhine from Peru (Fleagle and Tejedor 2002). It remains unclear how this older date was derived, however, and requires further substantiation. Atelidae | Ateles-Alouatta | 12.8 Ma–18 Ma | SAB The minimum age for the divergence of crown Atelidae (as recognized by Rosenberger 2011; subfamily Atelinae of others) is based on the first appearance of the crown atelid Stirtonia (Hershkovitz 1970; Rosenberger 2011) at 12.8 Ma. Fossil atelids are known from the middle Miocene of Colombia and the Quaternary of Brazil and the Greater Antilles (MacPhee and Horovitz 2002). The divergence of atelids is unlikely to have occurred before the first appearance of the potential stem or crown atelid Soriacebus at 18 Ma (Bown and Fleagle 1993; Fleagle et al. 1995). Cebidae | Cebus–Saimiri | 12.8 Ma–18 Ma | SAB The minimum age for the divergence of crown Cebidae is based on the first appearance of the crown cebid Neosaimiri (Stirton 1951; Hartwig and Meldrum 2002) dated at 12.8 from La Venta, Colombia. There are several older potential crown cebids including Dolichocebus and Tremacebus from Argentina and Chilecebus from Chile, all dated to around 20 Ma but it remains unclear how these taxa relate to the crown group. Recently, an additional potential crown cebid has been described from Panama (Bloch et al. 2016) dated at 20.9 Ma, which if substantiated would push the potential maximum bound to at least 21 Ma. Here we have accepted the notion that the divergence of cebids is unlikely to have occurred before the first appearance of the potential stem or crown atelid Soriacebus at 18 Ma (Bown and Fleagle 1993; Fleagle et al. 1995) but acknowledge that this date could be as old as 21–22 Ma. Cercopithecidae | Papionini–Cercopithecini | 5 Ma–23 Ma | SAB The minimum age for the divergence of crown Cercopithecidae is based on the first appearance of Parapapio in the late Miocene (5 Ma) at Lothagam, Kenya. It is potentially possible that some specimens of Parapapio from Lothagam could be as old as 7.4 Ma (Jablonski and Frost 2010). The divergence of cercopithecids is unlikely to have occurred before the first appearance of the stem cercopithecoid Prohylobates. The oldest documented Prohylobates specimens are from Wadi Moghra in the Qatarra Depression, Egypt dated to approximately 19.5 Ma. However, we view Kamoyapithecus, dated at 25 Ma as a crown hominoid, which implies that at least stem cercopithecoids were in existence at that time. Given the controversial nature of our views on Kamoyapithecus, we have used a date of 23 Ma as the likely maximum divergence time for crown cercopithecids. Colobinae | Colobini–Presbytini | 9.8 Ma–23 Ma | SAB The minimum age for the divergence of crown Colobinae is based on the first appearance of Microcolobus in the middle Miocene (9.8 Ma) at Ngeringerowa, Kenya (Benefit and Pickford 1986; Jablonski and Frost 2010). The divergence of colobines is unlikely to have occurred before the first appearance of the stem cercopithecoid Prohylobates. The oldest documented Prohylobates specimens are from Wadi Moghra in the Qatarra Depression, Egypt dated to approximately 19.5 Ma. However, we view Kamoyapithecus, dated at 25 Ma as a crown hominoid, which implies that at least stem cercopithecoids were in existence at that time. Given the controversial nature of our views on Kamoyapithecus, we have used a date of 23 Ma as the likely maximum divergence time for crown colobines. Appendix 2 Here, we briefly describe our new implementation of BF calculation in MCMCTree. Our approach is the thermodynamic integration-Gaussian quadrature method implemented recently in the program BPP. The mathematical details are given in Rannala and Yang (2017). In this article, we compare different rate models, $$m_{i}$$, which differ only in the density of the branch rates while the prior on divergence times and the sequence likelihood are the same between the models. The posterior distribution of times (t) and rates (r) given the sequence data $$D$$, and given a clock model $$m_{i}$$ is thus $f({\rm {\bf t}},{\rm {\bf r}}\vert D,m_{i} )=\frac{1}{z_{i} }f({\rm {\bf t}})f({\rm {\bf r}}\vert {\rm {\bf t}},m_{i} )f(D\vert {\rm {\bf t}},{\rm {\bf r}}),$ where $z_{i} =\int {f({\rm {\bf t}})f({\rm {\bf r}}\vert {\rm {\bf t}},m_{i} )f(D\vert {\rm {\bf t}},{\rm {\bf r}})\mbox{d}{\rm {\bf t}}\mbox{d}{\rm {\bf r}}}$ is the marginal likelihood of the data for model $$m_{i}$$. Let $$m_{r}$$ be the model with highest marginal likelihood, and let BF$$_{ir} \quad =$$$$z_{i}$$/$$z_{r}$$ be the BF of model $$i$$ over model $$r$$. Then the posterior probability of model $$i$$ is $\Pr (m_{i} \vert D)=\frac{z_{i} \Pr (m_{i} )}{\sum\limits_j {z_{j} \Pr (m_{j} )} }=\frac{z_{i} /z_{r} \Pr (m_{i} )}{\sum\limits_j {z_{j} /z_{r} \Pr (m_{j} )} }=\frac{\mbox{BF}_{ir} \Pr (m_{i} )}{\sum\limits_j {\mbox{BF}_{jr} \Pr (m_{j} )} },$ where the sum is over all models being tested, and Pr($$m_{\mathrm{i}})$$ is the prior model probability. We calculate $$z_{i}$$ by sampling from the power posterior $f_{\beta } ({\rm {\bf t}}\mbox{,}{\rm {\bf r}}\vert D,m_{i} )\propto f({\rm {\bf t}})f({\rm {\bf r}}\vert {\rm {\bf t}},m_{i} )f(D\vert {\rm {\bf t}}\mbox{,}{\rm {\bf r}})^{\beta },\quad 0\leqslant \beta \leqslant 1$ for a given $${\beta}$$ value. We choose $${K}\, \beta_{j}$$ values to integrate between 0 and 1 according to the Gauss-Legendre quadrature rule. The estimate of log $$z_i$$ is then given by the quadrature formula $\log z_{i} \approx \frac{1}{2}\sum\limits_{j=1}^K {w_{j} } \bar{{\ell }}_{\beta_{j} } ,$ where $$w_{j}$$ are the Gauss–Legendre quadrature weights, and $$\bar{{\ell }}_{\beta_{j} }$$ is the average of log-likelihood values, $$\ell_{\beta_{j}}$$, sampled from the power posterior with$$\beta_{j}$$. The standard error of the estimate is given by $\mbox{SE}=\frac{1}{2}\sqrt {\sum\limits_{j=1}^K {w_{j}^{2} \mbox{Var}} \left( {\ell_{\beta_{j} } } \right)/\mbox{ESS}} ,$ where the variance is calculated over the MCMC sample of $$\ell_{\beta_{j}}$$, and ESS is the effective-sample size of $$\ell_{\beta_{j}}$$. References Alfaro M.E., Santini F., Brock C., Alamillo H., Dornburg A., Rabosky D.L., Carnevale G., Harmon L.J. 2009 . Nine exceptional radiations plus high turnover explain species diversity in jawed vertebrates. Proc. Natl. Acad. Sci. USA 106 : 13410 – 13414 . Google Scholar CrossRef Search ADS Andersen M.J., Shult H.T., Cibois A., Thibault J.C., Filardi C.E., Moyle R.G. 2015 . Rapid diversification and secondary sympatry in Australo-Pacific kingfishers (Aves: Alcedinidae: Todiramphus). R. Soc. Open Sci., 2 : 140375 . Google Scholar CrossRef Search ADS Aris-Brosou S., Yang Z. 2002 . Effects of models of rate evolution on estimation of divergence dates with special reference to the metazoan 18S ribosomal RNA phylogeny. Syst. Biol. 51 : 703 – 714 . Google Scholar CrossRef Search ADS PubMed Badgley C., Behrensmeyer A.K. 1995 . Two long geological records of continental ecosystems. Palaeogeogr. Palaeoclimatol. Palaeoecol. 115 : 1 – 11 . Google Scholar CrossRef Search ADS Baele G., Lemey P., Bedford T., Rambaut A., Suchard M.A., Alekseyenko A.V. 2012a . Improving the accuracy of demographic and molecular clock model comparison while accommodating phylogenetic uncertainty. Mol. Biol. Evol. 29 : 2157 – 2167 . Google Scholar CrossRef Search ADS Baele G., Li W.L.S., Drummond A.J., Suchard M.A., Lemey P. 2012b . Accurate model selection of relaxed molecular clocks in Bayesian phylogenetics. Mol. Biol. Evol. 30 : 239 – 243 . Google Scholar CrossRef Search ADS Beard K. 1998 . A new genus of Tarsiidae (Mammalia: Primates) from the middle Eocene of Shanxi Province, China, with notes on the historical biogeography of tarsiers. Bull. Carnegie Mus. Nat. Hist. 34 : 260 – 277 . Beard K.C. 2016 . Out of Asia: Anthropoid origins and the colonization of Africa. Annu. Rev. Anthropol., 45 : 13.1 – 13.15 . Google Scholar CrossRef Search ADS Beard K.C., Qi T., Dawson M.R., Wang B., Li C. 1994 . A diverse new primate fauna from middle Eocene fissure-fillings in southeastern China. Nature 368 : 604 – 609 . Google Scholar CrossRef Search ADS PubMed Begun D.R. 2010 . Miocene hominids and the origins of the African apes and humans. Annu. Rev. Anthropol. 39 : 67 – 84 . Google Scholar CrossRef Search ADS Begun D.R., Nargolwalla M.C., Kordos L. 2012 . European Miocene hominids and the origin of the African ape and human clade. Evol. Anthropol. 21 : 10 – 23 . Google Scholar CrossRef Search ADS PubMed Benefit B.R., Pickford M. 1986 . Miocene fossil cercopithecoids from Kenya. Am. J. Phys. Anthropol. 69 : 441 – 464 . Google Scholar CrossRef Search ADS Benton M.J., Donoghue P.C. 2007 . Paleontological evidence to date the tree of life. Mol. Biol. Evol. 24 : 26 – 53 . Google Scholar CrossRef Search ADS PubMed Benton M.J., Donoghue P.C.J., Asher R.J. 2009 . Calibrating and constraining molecular clocks. In: Hedges B.S., Kumar S., editors. The timetree of life. Oxford, England : Oxford University Press . p. 35 – 86 . Bininda-Emonds O.R., Cardillo M., Jones K.E., MacPhee R.D., Beck R.M., Grenyer R., Price S.A., Vos R.A., Gittleman J.L., Purvis A. 2007 . The delayed rise of present-day mammals. Nature 446 : 507 – 512 . Google Scholar CrossRef Search ADS PubMed Bishop W. 1962 . The mammalian fauna and geomorphological relations of the Napak volcanics, Karamoja. Geol. Surv. Uganda Rec. 1957–1958 : 1 – 18 . Bloch J.I., Silcox M.T., Boyer D.M., Sargis E.J. 2007 . New Paleocene skeletons and the relationship of plesiadapiforms to crown-clade primates. Proc. Natl. Acad. Sci. USA 104 : 1159 – 1164 . Google Scholar CrossRef Search ADS Bloch J.I., Woodruff E.D., Wood A.R., Rincon A.F., Harrington A.R., Morgan, Gary S. Foster, David A. Camilo Montes, Jaramillo, Carlos A. Jud, Nathan A. Jones Douglas S. & MacFadden Bruce J. 2016 . First North American fossil monkey and early Miocene tropical biotic interchange. Nature 533 : 243 – 246 . Google Scholar CrossRef Search ADS PubMed Bond M., Tejedor M.F., Campbell K.E. Jr., Chornogubsky L., Novo N., Goin F. 2015 . Eocene primates of South America and the African origins of New World monkeys. Nature 520 : 538 – 541 . Google Scholar CrossRef Search ADS PubMed Bown T.M., Fleagle J.G. 1993 . Systematics, biostratigraphy, and dental evolution of the Palaeothentidae, later Oligocene to early-middle Miocene (Deseadan-Santacrucian) caenolestoid marsupials of South America. J. Paleontol. 2 Suppl : 1 – 76 . Brunet M., Guy F., Pilbeam D., Lieberman D.E., Likius A., Mackaye H.T., Ponce de Leon M.S., Zollikofer C.P., Vignaud P. 2005 . New material of the earliest hominid from the Upper Miocene of Chad. Nature 434 : 752 – 755 . Google Scholar CrossRef Search ADS PubMed Brunet M., Guy F., Pilbeam D., Mackaye H.T., Likius A., Ahounta D., Beauvilain A., Blondel C., Bocherens H., Boisserie J.-R., Louis De Bonis, Yves Coppens, Jean Dejax, Christiane Denys, Philippe Duringer, Véra Eisenmann, Gongdibé Fanone, Pierre Fronty, Denis Geraads, Thomas Lehmann, Fabrice Lihoreau, Antoine Louchart, Adoum Mahamat, Gildas Merceron, Guy Mouchelin, Olga Otero, Pablo Pelaez Campomanes, Marcia Ponce De Leon, Jean-Claude Rage, Michel Sapanet, Mathieu Schuster, Jean Sudre, Pascal Tassy, Xavier Valentin, Patrick Vignaud, Laurent Viriot, Antoine Zazzo, Christoph Zollikofer. 2002 . A new hominid from the upper Miocene of Chad, Central Africa. Nature 418 : 145 – 151 . Google Scholar CrossRef Search ADS PubMed Buckley G.A. 1997 . A new species of Purgatorius (Mammalia; Primatomorpha) from the lower Paleocene Bear formation, Crazy Mountains basin, south-central Montana. J. Paleontol. 71 : 149 – 155 . Google Scholar CrossRef Search ADS Burgess R., Yang Z. 2008 . Estimation of hominoid ancestral population sizes under Bayesian coalescent models incorporating mutation rate variation and sequencing errors. Mol. Biol. Evol. 25 : 1979 – 1994 . Google Scholar CrossRef Search ADS PubMed Chatterjee H.J., Ho S.Y., Barnes I., Groves C. 2009 . Estimating the phylogeny and divergence times of primates using a supermatrix approach. BMC Evol. Biol. 9 : 259 . Google Scholar CrossRef Search ADS PubMed Chester S.G., Bloch J.I., Boyer D.M., Clemens W.A. 2015 . Oldest known euarchontan tarsals and affinities of Paleocene Purgatorius to Primates. Proc. Natl. Acad. Sci. USA 112 : 1487 – 1492 . Clemens W.A. 1974 . Purgatorius, an early paromomyid primate (Mammalia). Science 184 : 903 – 905 . Google Scholar CrossRef Search ADS PubMed Clemens W.A. 2004 . Purgatorius (Plesiadapiformes, Primates, Mammalia), a Paleocene immigrant into northeastern Montana: stratigraphic occurrences and incisor proportions. Bull. Carnegie Mus. Nat. His. 36 : 3 – 13 . Google Scholar CrossRef Search ADS Corruccini R.S., Baba M., Goodman M., Ciochon R.L., Cronin J.E. 1980 . Non-linear macromolecular evolution and the molecular clock. Evolution 34 : 1216 – 1219 . Google Scholar CrossRef Search ADS PubMed de Queiroz A. 2005 . The resurrection of oceanic dispersal in historical biogeography. Trends Ecol. Evol. 20 : 68 – 73 Google Scholar CrossRef Search ADS PubMed de Queiroz A. 2014 . The monkey’s voyage: how improbable journeys shaped the history of life. New York : Basic Books . Dornburg A, Beaulieu J.M., Oliver J.C., Near T.J. 2011 . Integrating fossil preservation biases in the selection of calibrations for molecular divergence time estimation. Syst. Biol. 60 : 519 – 527 . Google Scholar CrossRef Search ADS PubMed dos Reis M., Donoghue P.C., Yang Z. 2014a . Neither phylogenomic nor palaeontological data support a Palaeogene origin of placental mammals. Biol. Lett. 10 : 20131003 . Google Scholar CrossRef Search ADS dos Reis M., Donoghue P.C., Yang Z. 2016 . Bayesian molecular clock dating of species divergences in the genomics era. Nat. Rev. Genet. 17 : 71 – 80 . Google Scholar CrossRef Search ADS PubMed dos Reis M., Inoue J., Hasegawa M., Asher R., Donoghue P.C., Yang Z. 2012 . Phylogenomic data sets provide both precision and accuracy in estimating the timescale of placental mammal phylogeny. Proc. R. Soc. Lond. B Biol. Sci. 279 : 3491 – 3500 . Google Scholar CrossRef Search ADS dos Reis M., Thawornwattana Y., Angelis K., Telford M.J., Donoghue P.C., Yang Z. 2015 . Uncertainty in the timing of origin of animals and the limits of precision in molecular timescales. Curr. Biol. 25 : 2939 – 2950 . Google Scholar CrossRef Search ADS PubMed dos Reis M., Yang Z. 2011 . Approximate likelihood calculation for Bayesian estimation of divergence times. Mol. Biol. Evol. 28 : 2161 – 2172 . Google Scholar CrossRef Search ADS PubMed dos Reis M., Yang Z. 2013 . The unbearable uncertainty of Bayesian divergence time estimation. J. Syst. Evol. 51 : 30 – 43 . Google Scholar CrossRef Search ADS dos Reis M., Zhu T., Yang Z. 2014b . The impact of the rate prior on Bayesian estimation of divergence times with multiple loci. Syst. Biol. 63 : 555 – 565 . Google Scholar CrossRef Search ADS Drummond A., Ho S., Phillips M., Rambaut A. 2006 . Relaxed phylogenetics and dating with confidence. PLoS Biol. 4 : e88 . Google Scholar CrossRef Search ADS PubMed Duchene S., Lanfear R., Ho S.Y.W. 2014 . The impact of calibration and clock-model choice on molecular estimates of divergence times. Mol. Phylogenet. Evol. 78 : 277 – 289 . Google Scholar CrossRef Search ADS PubMed Eizirik E., Murphy W.J., Springer M.S., O’Brien S.J. 2004 . Molecular phylogeny and dating of early primate divergences. In: Ross C.F., Kay R.F., editors. Anthropoid origins. Boston, MA : Springer. p. 45 – 64 . Google Scholar CrossRef Search ADS Finstermeier K., Zinner D., Brameier M., Meyer M., Kreuz E., Hofreiter M., Roos C. 2013 . A mitogenomic phylogeny of living primates. PLoS One 8 : e69504 . Google Scholar CrossRef Search ADS PubMed Fleagle J.G. 1990 . New fossil platyrrhines from the pinturas formation. South Argent. J. Hum. Evol. 19 : 61 – 85 . Google Scholar CrossRef Search ADS Fleagle J.G., Bown T., Swisher C., Buckley G. 1995 . Age of the Pinturas and Santa Cruz formations. VI Congreso Argentino de Paleontología y Bioestratigrafía. Trelew, Chubut, Argentine : Museo Paleontológico Eugidio Feruglio p. 129 – 135 . Fleagle J.G., Powers D.W., Conroy G.C., Watters J.P. 1987 . New fossil platyrrhines from Santa-Cruz Province, Argentina. Folia Primatol. 48 : 65 – 77 . Google Scholar CrossRef Search ADS Fleagle J.G., Tejedor M.F. 2002 . “Early platyrrhines of southern South America” in Hartwig, W.C. (Ed). The Primate fossil record. Cambridge : Cambridge University Press, pp. 161 – 174 . Flynn L.J., Barry J.C., Morgan M.E., Pilbeam D., Jacobs L.L., Lindsay E.H. 1995 . Neogene Siwalik mammalian lineages—species longevities, rates of change, and modes of speciation. Palaeogeogr. Palaeocl. 115 : 249 – 264 . Google Scholar CrossRef Search ADS Fox R.C., Scott C.S. 2011 . A new, early Puercan (Earliest Paleocene) species of Purgatorius (Plesiadapiformes, Primates) from Saskatchewan, Canada. J. Paleontol. 85 : 537 – 548 . Google Scholar CrossRef Search ADS Gelman A., Carlin J.B., Stern H.S., Dunson D.B., Vehtari A., Rubin D.B. 2013 . Bayesian data analysis, Vol. 21, 3rd ed. Boca Raton, FL : CRC Press . Geraads D., Alemseged Z., Bellon H. 2002 . The late Miocene mammalian fauna of Chorora, Awash basin, Ethiopia: systematics, biochronology and 40K-40Ar ages of the associated volcanics. Tertiary Res. 113 – 122 . Gheerbrant E., Cappetta H., Feist M., Jaeger J-J., Sudre J., Vianey-Liaud M., Sigé B. 1993 . La succession des faunes de vertébrés d’âge paléocène supérieur et éocène inférieur dans le bassin d’Ouarzazate, Maroc. Contexte géologique, portée biostratigraphique et paléogéographique. Newsl Stratigr 33 – 58 . Gheerbrant E., Sudre J., Sen S., Abrial C., Marandat B., Sigé B., Vianey-Liaud M. 1998 . Nouvelles données sur les mammifères du Thanétien et de l’Yprésien du Bassin d’Ouarzazate (Maroc) et leur contexte stratigraphique. Palaeovertebrata 27 : 155 – 202 . Gingerich P.D. 1981 . Eocene Adapidae, paleobiogeography, and the origin of South American Platyrrhini. In: Ciochon R.L., Chiarelli A.B., editors. Evolutionary biology of the new world monkeys and continental drift. New York : Plenum Publishing Corporation. pp. 123 – 138 . Google Scholar CrossRef Search ADS Goodman M., Braunitzer G., Stangl A., Schrank B. 1983 . Evidence on human origins from haemoglobins of African apes. Nature 303 : 546 – 548 . Google Scholar CrossRef Search ADS PubMed Goodman M., Porter C.A., Czelusniak J., Page S.L., Schneider H., Shoshani J., Gunnell G., Groves C.P. 1998 . Toward a phylogenetic classification of Primates based on DNA evidence complemented by fossil evidence. Mol. Phylogenet. Evol. 9 : 585 – 598 . Google Scholar CrossRef Search ADS PubMed Graur D., Martin W. 2004 . Reading the entrails of chickens: molecular timescales of evolution and the illusion of precision. Trends Genet. 20 : 80 – 86 . Google Scholar CrossRef Search ADS PubMed Harrison T. 2010a . Dendropithecoidea, proconsuloidea, and hominoidea. In: Werdelin L., Sanders W., editors. Cenozoic mammals of Africa. Berkeley, CA : University of California Press. pp. 429 – 469 . Google Scholar CrossRef Search ADS Harrison T. 2010b . Later Tertiary lorisiformes. In: Werdelin L., Sanders W., editors. Cenozoic mammals of Africa. Berkeley, CA : University of California Press . pp. 333 – 349 . Google Scholar CrossRef Search ADS Hartwig W.C., Meldrum D.J. 2002 . Miocene platyrrhines of the northern Neotropics. Cam. S. Biol. Evol. Anthropol. 175 – 188 . Hasegawa M., Kishino H., Yano T. 1985 . Dating the human-ape splitting by a molecular clock of mitochondrial DNA. J. Mol. Evol. 22 : 160 – 174 . Google Scholar CrossRef Search ADS PubMed Heath T.A., Holder M.T., Huelsenbeck J.P. 2012 . A Dirichlet process prior for estimating lineage-specific substitution rates. Mol. Biol. Evol. 29 : 939 – 955 . Google Scholar CrossRef Search ADS PubMed Heath T.A., Huelsenbeck J.P., Stadler T. 2014 . The fossilized birth–death process for coherent calibration of divergence-time estimates. Proc. Natl. Acad. Sci. USA 111 : E2957 – E2966 . Google Scholar CrossRef Search ADS Hershkovitz P. 1970 . Notes on tertiary platyrrhine monkeys and description of a new genus from the late Miocene of Colombia. Folia Primatol. 12 : 1 – 37 . Google Scholar CrossRef Search ADS PubMed Ho S.Y. 2009 . An examination of phylogenetic models of substitution rate variation among lineages. Biol. Lett. 5 : 421 – 424 . Google Scholar CrossRef Search ADS PubMed Ho S.Y., Duchene S. 2014 . Molecular-clock methods for estimating evolutionary rates and timescales. Mol. Ecol., 23 : 5947 – 5965 . Google Scholar CrossRef Search ADS Ho S.Y., Duchene S., Duchene D. 2015 . Simulating and detecting autocorrelation of molecular evolutionary rates among lineages. Mol. Ecol. Resour. 15 : 688 – 696 . Google Scholar CrossRef Search ADS PubMed Ho S.Y., Philips M.J. 2009 . Accounting for calibration uncertainty in phylogenetic estimation of evolutionary divergence times. Syst. Biol. 58 : 367 – 380 . Google Scholar CrossRef Search ADS PubMed Ho S.Y., Phillips M.J., Drummond A.J., Cooper A. 2005 . Accuracy of rate estimation using relaxed-clock models with a critical focus on the early metazoan radiation. Mol. Biol. Evol. 22 : 1355 – 1363 . Google Scholar CrossRef Search ADS PubMed Ho, S.Y., 2014 . The changing face of the molecular clock. Trends Ecol Evol. 29 : 496 – 503 . Google Scholar CrossRef Search ADS PubMed Inoue J., Donoghue P.C.H., Yang Z. 2010 . The impact of the representation of fossil calibrations on Bayesian estimation of species divergence times. Syst. Biol. 59 : 74 – 89 . Google Scholar CrossRef Search ADS PubMed Jablonski N., Frost S. 2010 . Cercopithecoidea. In: Werdelin L., Sanders. W., editors. Cenozoic mammals of Africa. Berkeley : University of California Press . p. 393 – 428 . Google Scholar CrossRef Search ADS Jacobs L.L. 1981 . Miocene lorisid primates from the Pakistan Siwaliks. Nature 289 : 585 – 587 . Google Scholar CrossRef Search ADS PubMed Kappelman J., Kelley J., Pilbeam D., Sheikh K.A., Ward S., Anwar M., Barry J.C., Brown B., Hake P., Johnson N.M. et al. 1991 . The earliest occurrence of Sivapithecus from the middle Miocene Chinji Formation of Pakistan. J. Hum. Evol. 21 : 61 – 73 . Google Scholar CrossRef Search ADS Kay R.F., Johnson D., Meldrum D.J. 1998 . A new pitheciin primate from the middle Miocene of Argentina. Am. J. Primatol., 45 : 317 – 336 . Google Scholar CrossRef Search ADS Kishino H., Thorne J.L., Bruno W.J. 2001 . Performance of a divergence time estimation method under a probabilistic model of rate evolution. Mol. Biol. Evol. 18 : 352 – 361 . Google Scholar CrossRef Search ADS PubMed Koop B.F., Tagle D.A., Goodman M., Slightom J.L. 1989 . A molecular view of primate phylogeny and important systematic and evolutionary questions. Mol. Biol. Evol. 6 : 580 – 612 . Google Scholar PubMed Kumar S., Hedges S.B. 1998 . A molecular timescale for vertebrate evolution. Nature 392 : 917 – 920 . Google Scholar CrossRef Search ADS PubMed Lartillot N., Philippe H. 2006 . Computing Bayes factors using thermodynamic integration. Syst. Biol. 55 : 195 – 207 . Google Scholar CrossRef Search ADS PubMed Lartillot N., Phillips M.J., Ronquist F. 2016 . A mixed relaxed clock model. Phil. Trans. R. Soc. Lond. B 371 : 20150132 . Google Scholar CrossRef Search ADS Leakey M.G., Ungar P.S., Walker A. 1995 . A new genus of large primate from the late Oligocene of Lothidok, Turkana District, Kenya. J. Hum. Evol. 28 : 519 – 531 . Google Scholar CrossRef Search ADS Lebatard A.E., Bourles D.L., Duringer P., Jolivet M., Braucher R., Carcaillet J., Schuster M., Arnaud N., Monie P., Lihoreau F. et al. 2008 . Cosmogenic nuclide dating of Sahelanthropus tchadensis and Australopithecus bahrelghazali: Mio-Pliocene hominids from Chad. Proc. Natl. Acad. Sci. USA 105 : 3226 – 3231 . Google Scholar CrossRef Search ADS Lepage T., Bryant D., Philippe H., Lartillot N. 2007 . A general comparison of relaxed molecular clock models. Mol. Biol. Evol. 24 : 2669 – 2680 . Google Scholar CrossRef Search ADS PubMed Linder H.P., Hardy C.R., Rutschmann F. 2005 . Taxon sampling effects in molecular clock dating: an example from the African Restionaceae. Mol. Phylogenet. Evol. 35 : 569 – 582 . Google Scholar CrossRef Search ADS PubMed Linder M., Britton T., Sennblad B. 2011 . Evaluation of Bayesian models of substitution rate evolution: parental guidance versus mutual independence. Syst. Biol. 60 : 329 – 342 . Google Scholar CrossRef Search ADS PubMed Liu Y., Medina R., Goffinet B. 2014 . 350 my of mitochondrial genome stasis in mosses, an early land plant lineage. Mol. Biol. Evol. 31 : 2586 – 2591 . Google Scholar CrossRef Search ADS PubMed Lovejoy C.O., Burstein H., Heiple K.H. 1972 . Primate phylogeny and immunological distance. Science 176 : 803 – 805 . Google Scholar CrossRef Search ADS PubMed Luo Z.X. 2007 . Transformation and diversification in early mammal evolution. Nature 450 : 1011 – 1019 . Google Scholar CrossRef Search ADS PubMed Luo Z.X., Yuan C.X., Meng Q.J., Ji Q. 2011 . A Jurassic eutherian mammal and divergence of marsupials and placentals. Nature 476 : 442 – 445 . Google Scholar CrossRef Search ADS PubMed McDougall I.A.N., Feibel C.S. 2003 . Numerical age control for the Miocene-Pliocene succession at Lothagam, a hominoid-bearing sequence in the northern Kenya Rift. In: Leakey M.G.H., Harris J.M., editors. Lothagam: the dawn of humanity in Eastern Africa. New York : Columbia University Press. p. 43 – 64 . MacInnes D.G. 1943 . Notes on the East African primates. J. East Afr. Uganda Nat. Hist. Soc. 39 : 521 – 530 . MacLatchy L., DeSilva J., Sanders W., Wood B. 2010 . Hominini. In: Werdelin L., Sanders W. editors. Cenozoic mammals of Africa. Berkeley : University of California Press . p. 471 – 540 . Google Scholar CrossRef Search ADS MacPhee R.D., Horovitz I. 2002 . Extinct quaternary platyrrhines of the Greater Antilles and Brazil. Cam. S. Biol. Evol. Anthropol. 189 – 200 . Madden C.T. 1980 . NewProconsul (Xenopithecus) from the Miocene of Kenya. Primates 21 : 241 – 252 . Google Scholar CrossRef Search ADS Marshall C.R. 1990 . The fossil record and estimating divergence times between lineages: maximum divergence times and the importance of reliable phylogenies. J. Mol. Evol. 30 : 400 – 408 . Google Scholar CrossRef Search ADS PubMed Marshall C.R. 2008 . A simple method for bracketing absolute divergence times on molecular phylogenies using multiple fossil calibration points. Am. Nat. 171 : 726 – 742 . Google Scholar CrossRef Search ADS PubMed Martin R.D. 1993 . Primate origins: plugging the gaps. Nature 363 : 223 – 234 . Google Scholar CrossRef Search ADS PubMed Martin R.D., Soligo C., Tavaré S. 2007 . Primate origins: implications of a Cretaceous ancestry. Folia Primatol. 78 : 277 – 296 . Google Scholar CrossRef Search ADS PubMed McCrossin M.L. 1999 . Phylogenetic relationships and paleoecological adaptations of a new bushbaby from the middle Miocene of Kenya. Am. J. Phys. Anthropol. 89 : 195 – 196 . Nascimento F.F., Lazar A., SeuÃ¡nez H.N., Bovincino C.R. 2014 . Reanalysis of the biogeographical hypothesis of range expansion between robust and gracile capuchin monkeys. J. Phylogeogr. 42 : 1349 – 1357 . Nascimento F.F., dos Reis M., Yang Z. 2017 . A biologist’s guide to Bayesian phylogenetic analysis. Nat. Ecol. Evol. 1 : 1446 – 1454 . Google Scholar CrossRef Search ADS PubMed O’Leary M.A., Bloch J.I., Flynn J.J., Gaudin T.J., Giallombardo A., Norberto P. Giannini, Suzann L. Goldberg, Kraatz, Brian P. Zhe-Xi Luo, Jin Meng, Xijun Ni, Novacek, Michael J. Perini, Fernando A. Randall, Zachary S. Rougier, Guillermo W. Sargis, Eric J. Silcox, Mary T. Simmons, Nancy B. Michelle Spaulding, Velazco, Paúl M. Marcelo Weksler, Wible, John R. Cirranello. Andrea L. 2013 . The placental mammal ancestor and the post-K-Pg radiation of placentals. Science 339 : 662 – 667 . Google Scholar CrossRef Search ADS PubMed O’Reilly J., dos Reis M., Donoghue P.C.J. 2015 . Dating tips for divergence time estimation. Trends Genet. 31 : 637 – 650 . Google Scholar CrossRef Search ADS PubMed Paradis E. 2013 . Molecular dating of phylogenies by likelihood methods: a comparison of models and a new information criterion. Mol. Phylogenet. Evol. 67 : 436 – 444 . Google Scholar CrossRef Search ADS PubMed Perelman P., Johnson W.E., Roos C., SeuÃ¡nez H.N., Horvath J.E., Moreira M.A.M., Kessing B., Pontius J., Roelke M., Rumpler Y., Schneider, Maria Paula C. Artur Silva, O’Brien, Stephen J. Pecon-Slattery. Jill 2011 . A molecular phylogeny of living primates. PLoS Genet. 7 : e1001342 EP . Google Scholar CrossRef Search ADS PubMed Phillips E.M., Walker A. 2000 . A new species of fossil lorisid from the miocene of east Africa. Primates 41 : 367 – 372 . Google Scholar CrossRef Search ADS Phillips E.M., Walker A. 2002 . Fossil lorisoids. In: Hartwig W.C., editor. The primate fossil record. Cambridge : Cambridge University Press. p. 83 – 95 . Poux C., Douzery E.J. 2004 . Primate phylogeny, evolutionary rate variations, and divergence times: a contribution from the nuclear gene IRBP. Am. J. Phys. Anthropol. 124 : 1 – 16 . Google Scholar CrossRef Search ADS PubMed Pozzi L., Hodgson J.A., Burrell A.S., Sterner K.N., Raaum R.L., Disotell T.R. 2014 . Primate phylogenetic relationships and divergence dates inferred from complete mitochondrial genomes. Mol. Phylogenet. Evol. 75 : 165 – 183 . Google Scholar CrossRef Search ADS PubMed Prum R.O., Berv J.S., Dornburg A., Field D.J., Townsend J.P., Lemmon E.M., Lemmon A.R. 2015 . A comprehensive phylogeny of birds (Aves) using targeted next-generation DNA sequencing. Nature 526 : 569 – 573 . Google Scholar CrossRef Search ADS PubMed Purvis A. 1995 . A composite estimate of primate phylogeny. Phil. Trans. R. Soc. Lond. B , 348 : 405 – 421 . Google Scholar CrossRef Search ADS Radinsky L. 1978 . Do albumin clocks run on time? Science 200 : 1182 – 1183 . Google Scholar CrossRef Search ADS PubMed Rannala B. 2016 . Conceptual issues in Bayesian divergence time estimation. Phil. Trans. R. Soc. B, 371 : 20150134 . Google Scholar CrossRef Search ADS Rannala B., Yang Z. 2007 . Inferring speciation times under an episodic molecular clock. Syst. Biol. 56 : 453 – 466 . Google Scholar CrossRef Search ADS PubMed Rannala B., Yang Z. 2017 . Efficient Bayesian species tree inference under the multispecies coalescent. Syst. Biol. 66 : 823 – 842 . Google Scholar PubMed Rasmussen D.T., Gutierrez M. 2009 . A mammalian fauna from the late oligocene of Northwestern Kenya. Palaeontogr. Abt. A 288 : 1 – 52 . Raza S.M., Barry J.C., Pilbeam D., Rose M.D., Shah S.M.I., Ward S. 1983 . New hominoid primates from the middle Miocene Chinji Formation, Potwar Plateau, Pakistan. Nature 306 : 52 – 54 . Google Scholar CrossRef Search ADS Read D.W., Lestrel P.E. 1970 . Hominid phylogeny and immunology: a critical appraisal. Science 168 : 578 – 580 . Google Scholar CrossRef Search ADS PubMed Rieux A., Eriksson A., Li M., Sobkowiak B., Weinert L.A., Warmuth V., Ruiz-Linares A., Manica A., Balloux F. 2014 . Improved calibration of the human mitochondrial clock using ancient genomes. Mol. Biol. Evol. 31 : 2780 – 2792 . Google Scholar CrossRef Search ADS PubMed Ronquist F., Klopfstein S., Vilhelmsen L., Schulmeister S., Murray D.L., Rasnitsyn P. 2012 . A total-evidence approach to dating with fossils, applied to the early radiation of the Hymenoptera. Syst. Biol. 61 : 973 – 999 . Google Scholar CrossRef Search ADS PubMed Rosenberger A. 2011 . Evolutionary morphology, platyrrhine evolution, and systematics. Anat. Rec. 294 : 1955 – 1974 . Google Scholar CrossRef Search ADS Rutschmann F., Eriksson T., Abu Salim K., Conti E. 2007 . Assessing calibration uncertainty in molecular dating: the assignment of fossils to alternative calibration points. Syst. Biol. 56 : 591 – 608 . Google Scholar CrossRef Search ADS PubMed Sanderson M.J. 1997 . A nonparametric approach to estimating divergence times in the absence of rate constancy. Mol. Biol. Evol. 14 : 1218 – 1231 . Google Scholar CrossRef Search ADS Sarich V.M., Wilson A.C. 1967 . Immunological time scale for Hominoid evolution. Science 158 : 1200 – 1203 Google Scholar CrossRef Search ADS PubMed Seiffert E.R. 2010 . Chronology of Paleogene mammal localities. In: Werdelin L., Sanders W., editors. Cenozoic mammals of Africa. Berkeley : University of California Press. p. 19 – 26 . Google Scholar CrossRef Search ADS Seiffert E.R. 2012 . Early primate evolution in Afro-Arabia. Evol. Anthropol. 21 : 239 – 253 . Google Scholar CrossRef Search ADS PubMed Seiffert E.R., Simons E.L., Attia Y. 2003 . Fossil evidence for an ancient divergence of lorises and galagos. Nature 422 : 421 – 424 . Google Scholar CrossRef Search ADS PubMed Shapiro B., Rambaut A., Drummond A.J. 2006 . Choosing appropriate substitution models for the phylogenetic analysis of protein-coding sequences. Mol. Biol. Evol. 23 : 7 – 9 . Google Scholar CrossRef Search ADS PubMed Simpson G.G. 1967 . The tertiary lorisiform primates of Africa. Bull. Mus. Comp. Zool. 136 : 39 – 62 . Silcox M.T. 2008 . The biogeographic origins of primates and euprimates: East, West, North, or South of Eden? In: Dagosto M., Sargis E.J., editors. Mammalian evolutionary morphology: a tribute to Frederick S. Szalay. New York : Springer . p. 199 – 231 . Google Scholar CrossRef Search ADS Soltis P.S., Soltis D.E., Savolainen V., Crane P.R., Barraclough T.G. 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 . Google Scholar CrossRef Search ADS Springer M.S., Meredith R.W., Gatesy J., Emerling C.A., Park J., Rabosky D.L., Stadler T., Steiner C., Ryder O.A, Janečka J.E, Fisher, Colleen A. Murphy. William J. 2012 . Macroevolutionary dynamics and historical biogeography of primate diversification inferred from a species supermatrix. PLoS One 7 : e49521 . Google Scholar CrossRef Search ADS PubMed Stamatakis A. 2014 . RAxML version 8: a tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics 30 : 1312 – 1313 . Google Scholar CrossRef Search ADS PubMed Steiper M.E., Young N.M. 2006 . Primate molecular divergence dates. Mol. Phylogenet. Evol. 41 : 384 – 394 . Google Scholar CrossRef Search ADS PubMed Steiper M.E., Seiffert E.R. 2012 . Evidence for a convergent slowdown in primate molecular rates and its implications for the timing of early primate evolution. Proc. Nat. Acad. Sci. USA 109 : 6006 – 6011 . Google Scholar CrossRef Search ADS Stirton R.A. 1951 . Ceboid monkeys from the Miocene of Colombia. Bull. Univ. Calif. Pub. Geol. Sci. 28 : 315 – 356 . Suwa G., Kono R.T., Katoh S., Asfaw B., Beyene Y. 2007 . A new species of great ape from the late Miocene epoch in Ethiopia. Nature 448 : 921 – 924 . Google Scholar CrossRef Search ADS PubMed Thorne J.L., Kishino H. 2002 . Divergence time and evolutionary rate estimation with multilocus data. Syst. Biol. 51 : 689 – 702 . Google Scholar CrossRef Search ADS PubMed Thorne J.L., Kishino H., Painter I.S. 1998 . Estimating the rate of evolution of the rate of molecular evolution. Mol. Biol. Evol. 15 : 1647 – 1657 . Google Scholar CrossRef Search ADS PubMed Uzzell T., Pilbeam D. 1971 . Phyletic divergence dates of hominoid primates: a comparison of fossil and molecular data. Evolution 25 : 615 – 635 . Google Scholar CrossRef Search ADS PubMed Walker A.C. 1978 . Prosimian primates. In: Maglio V.J., Cooke H.B.S., editors. Evolution of African mammals. Cambridge : Harvard University Press. p. 90 – 99 . Warnock R.C., Parham J.F., Joyce W.G., Lyson T.R., Donoghue P.C. 2015 . Calibration uncertainty in molecular dating analyses: there is no substitute for the prior evaluation of time priors. Proc. R. Soc. B 282 : 20141013 . Google Scholar CrossRef Search ADS Wilkinson R.D., Steiper M.E., Soligo C., Martin R.D., Yang Z., Tavare S. 2011 . Dating primate divergences through an integrated analysis of palaeontological and molecular data. Syst. Biol. 60 : 16 – 31 . Google Scholar CrossRef Search ADS PubMed Wolpoff M.H., Hawks J., Senut B., Pickford M., Ahern J. 2006 . An ape or the ape: is the Toumaï cranium TM 266 a hominid. Paleoanthropology 2006 : 36 – 50 . Xie W., Lewis P.O., Fan Y., Kuo L., Chen M.-H. 2011 . Improving marginal likelihood estimation for Bayesian phylogenetic model selection. Syst. Biol. 60 : 150 – 160 . Google Scholar CrossRef Search ADS PubMed Yang Z. 1994a . Estimating the pattern of nucleotide substitution. J. Mol. Evol. 39 : 105 – 111 . Yang Z. 1994b . Maximum likelihood phylogenetic estimation from DNA sequences with variable rates over sites: approximate methods. J. Mol. Evol. 39 : 306 – 314 . Google Scholar CrossRef Search ADS Yang Z. 2007 . PAML 4: phylogenetic analysis by maximum likelihood. Mol. Biol. Evol. 24 : 1586 – 1591 . Google Scholar CrossRef Search ADS PubMed Yang Z., Rannala B. 1997 . Bayesian phylogenetic inference using DNA sequences: a Markov chain Monte Carlo Method. Mol. Biol. Evol. 14 : 717 – 724 . Google Scholar CrossRef Search ADS PubMed Yang Z., Rannala B. 2006 . Bayesian estimation of species divergence times under a molecular clock using multiple fossil calibrations with soft bounds. Mol. Biol. Evol. 23 : 212 – 226 . Google Scholar CrossRef Search ADS PubMed Yang Z., Yoder A.D. 2003 . Comparison of likelihood and Bayesian methods for estimating divergence times using multiple gene loci and calibration points, with application to a radiation of cute-looking mouse lemur species. Syst. Biol. 52 : 705 – 716 . Google Scholar CrossRef Search ADS PubMed Yoder A.D., Nowak M.D. 2006 . Has vicariance or dispersal been the predominant biogeographic force in Madagascar? Only time will tell. Annu. Rev. Ecol. Evol. Syst. 37 : 404 – 431 . Google Scholar CrossRef Search ADS Yoder A.D., Yang Z. 2004 . Divergence dates for Malagasy lemurs estimated from multiple gene loci: geological and evolutionary context. Mol. Ecol. 13 : 757 – 773 . Google Scholar CrossRef Search ADS PubMed Zalmout I.S., Sanders W.J., Maclatchy L.M., Gunnell G.F., Al-Mufarreh Y.A., Ali M.A., Nasser A.A., Al-Masari A.M., Al-Sobhi S.A., Nadhra A.O, Matari, Adel H. Wilson Jeffrey A. & Gingerich. Philip D. 2010 . New Oligocene primate from Saudi Arabia and the divergence of apes and Old World monkeys. Nature 466 : 360 – 364 . Google Scholar CrossRef Search ADS PubMed Zhang C., Stadler T., Klopfstein S., Heath T.A., Ronquist F. 2016 . Total-evidence dating under the fossilized birth-death process. Syst. Biol. 65 : 228 – 249 . Google Scholar CrossRef Search ADS PubMed Zhou X., Xu S., Xu J., Chen B., Zhou K., Yang G. 2012 . Phylogenomic analysis resolves the interordinal relationships and rapid diversification of the laurasiatherian mammals. Syst. Biol. 61 : 150 – 164 . Google Scholar CrossRef Search ADS PubMed Zhu T, dos Reis M, Yang Z. 2015 . Characterization of the uncertainty of divergence time estimation under relaxed molecular clock models using multiple loci. Syst Biol 64 : 267 – 280 . Google Scholar CrossRef Search ADS PubMed Zuckerkandl E., Pauling L. 1965 . Molecules as documents of evolutionary history. J. Theor. Biol. 8 : 357 – 366 . Google Scholar CrossRef Search ADS PubMed © The Author(s) 2018. Published by Oxford University Press, on behalf of the Society of Systematic Biologists. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited. For Permissions, please email: journals.permissions@oup.com http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Systematic Biology Oxford University Press

# Using Phylogenomic Data to Explore the Effects of Relaxed Clocks and Calibration Strategies on Divergence Time Estimation: Primates as a Test Case

, Volume Advance Article (4) – Jan 12, 2018
22 pages

/lp/ou_press/using-phylogenomic-data-to-explore-the-effects-of-relaxed-clocks-and-c50FkNfg5A
Publisher
Oxford University Press
© The Author(s) 2018. Published by Oxford University Press, on behalf of the Society of Systematic Biologists.
ISSN
1063-5157
eISSN
1076-836X
D.O.I.
10.1093/sysbio/syy001
Publisher site
See Article on Publisher Site

### Abstract

Abstract Primates have long been a test case for the development of phylogenetic methods for divergence time estimation. Despite a large number of studies, however, the timing of origination of crown Primates relative to the Cretaceous–Paleogene (K–Pg) boundary and the timing of diversification of the main crown groups remain controversial. Here, we analysed a data set of 372 taxa (367 Primates and 5 outgroups, 3.4 million aligned base pairs) that includes nine primate genomes. We systematically explore the effect of different interpretations of fossil calibrations and molecular clock models on primate divergence time estimates. We find that even small differences in the construction of fossil calibrations can have a noticeable impact on estimated divergence times, especially for the oldest nodes in the tree. Notably, choice of molecular rate model (autocorrelated or independently distributed rates) has an especially strong effect on estimated times, with the independent rates model producing considerably more ancient age estimates for the deeper nodes in the phylogeny. We implement thermodynamic integration, combined with Gaussian quadrature, in the program MCMCTree, and use it to calculate Bayes factors for clock models. Bayesian model selection indicates that the autocorrelated rates model fits the primate data substantially better, and we conclude that time estimates under this model should be preferred. We show that for eight core nodes in the phylogeny, uncertainty in time estimates is close to the theoretical limit imposed by fossil uncertainties. Thus, these estimates are unlikely to be improved by collecting additional molecular sequence data. All analyses place the origin of Primates close to the K–Pg boundary, either in the Cretaceous or straddling the boundary into the Palaeogene. Bayes factors, Bayesian analysis, fossil, molecular dating, phylogenomic analysis, Primates, relaxed clock Divergence time estimation is fundamentally important to every field of evolutionary analysis. Reliable estimates of the timing of speciation events across the Tree of Life allow us to correlate these events with both biotic and abiotic phenomena on geological and more recent timescales, thus illuminating those that are most closely associated with periods of rapid diversification (Zhou et al. 2012; Andersen et al. 2015), evolutionary stasis (Alfaro et al. 2009; Liu et al. 2014), or high rates of extinction (Prum et al. 2015; Zhang et al. 2016). The prospect of employing genomic information for discovering the geological age of the Primates began virtually simultaneously with the publication of Zuckerkandl and Pauling’s (1965) molecular clock hypothesis. Sarich and Wilson (1967) employed a strict clock interpretation of immunological distance data to hypothesize that humans and other African apes (i.e., chimp and gorilla) shared a common ancestor as recently as 5 million years ago. This was revolutionary at the time given the implications for the necessarily rapid evolution of bipedal locomotion in the hominin lineage, and accordingly, drew considerable attention from the anthropological community (Read and Lestrel 1970; Uzzell and Pilbeam 1971; Lovejoy et al. 1972; Radinsky 1978; Corruccini et al. 1980). Despite this interest, it wasn’t until the 1980s that the field of divergence time estimation assumed a relatively modern flavor. It was then that investigators began to apply statistical models to DNA sequence data for the purposes of branch length and divergence time estimation (e.g., Hasegawa et al. 1985). Remarkably, these studies first emerged at a time when the sister–lineage relationship of humans to chimps was considered highly controversial (e.g., Goodman et al. 1983)—a relationship that is now considered unequivocal. Though it was Zuckerkandl and Pauling’s (1965) “molecular clock” hypothesis that made divergence time estimation feasible, it became clear very soon thereafter that there are myriad violations to a uniform clock. Thus, in subsequent decades, increasingly sophisticated models have been developed for relaxing the assumptions of a strict clock (Thorne et al. 1998; Kishino et al. 2001; Thorne and Kishino 2002; Drummond et al. 2006; Lepage et al. 2007; Rannala and Yang 2007; Heath et al. 2012). These models can be loosely divided into two categories: autocorrelated models, wherein rates of evolution in daughter species are statistically distributed around the parental rates (Sanderson 1997; Thorne et al. 1998; Aris-Brosou and Yang 2002), and uncorrelated models, wherein each lineage on the tree is free to assume a fully independent rate (Drummond et al. 2006; Rannala and Yang 2007; Paradis 2013). A parallel challenge to divergence time analysis can be observed in the development of calibration strategies (Marshall 1990; Yang and Rannala 2006; Benton and Donoghue 2007; Marshall 2008; Dornburg et al. 2011). Given that branch lengths on a phylogeny are the product of rate and time, investigators must assume one to infer the other. The most typical method for breaking this impasse is to employ fossil data as calibrations for one or more nodes in a given phylogeny so that the ages of other nodes can be inferred using DNA sequence data. This places an enormous burden on both the correct placement and age assignment of the fossils. If they are misplaced (i.e., assigned to clades where they do not belong) or if their geological ages are misinterpreted, divergence time estimates for the entire tree can be severely compromised (Martin 1993). The uncertainty imposed by paleontological ambiguity has not been as widely appreciated as have been the uncertainties introduced by the finite amount of DNA sequence data, which with the “genomics revolution,” is becoming steadily less problematic. We have reached a state of the art where branch lengths can be estimated with very high precision. The combination of genome-scale data, sophisticated molecular evolutionary models, and ever-increasing computational power has brought us to this point. Advances in DNA sequencing technology are such that virtually every major clade has at least a few species represented by whole-genome sequences, and this trend is rapidly accelerating. Bayesian methods have been developed such that parameter space can be explored during MCMC estimation, and though violations of the molecular clock will continue to present problems, methods for measuring and accommodating rate variation across phylogenies are explicit and generalizable. And finally, the computing power to integrate this information is increasing steadily. But because of the confounding effect of nonindependence of rate and time, the expectation of a conventional Bayesian analysis—that infinite data will eventually overcome prior information and provide posterior distributions with certainty—cannot be met (dos Reis and Yang 2013; Zhu et al. 2015). Thus, the field at present is focused on developing a better understanding of the effects of relaxed clock model choice, and on the impacts of calibration points, both with regard to abundance and placement in the phylogeny. Furthermore, in addressing these challenges, it is an open question as to whether simulation studies or tests of empirical data will be more informative for our understanding of best practices. With regard to clock model choice, an empirical study of three independent data sets showed that autocorrelated models outperform uncorrelated models, though the same study found a “high sensitivity” to the relaxation model employed (Lepage et al. 2007), while another empirical study found, however, that an independent rates model was superior (Linder et al. 2011). Simulation studies have only recently been employed, finding that even with complete taxon sampling, rate autocorrelation is challenging to detect (Ho et al. 2015). This has led to the conclusion that rigorous model selection should be conducted on a case-by-case basis, utilizing a wide range of real data sets, and thus comprising a promising avenue for future research (Ho et al. 2005; Ho and Duchene 2014; Ho et al. 2015). With regard to fossil calibration strategies, simulation studies (e.g., Duchene et al. 2014) have thus far agreed with previous empirical studies in finding that multiple calibrations are fundamentally important to accurate age estimation (Soltis et al. 2002; Yoder and Yang 2004; Linder et al. 2005; Rutschmann et al. 2007; Marshall 2008). Duchene et al. (2014) noted that calibrations close to the root of the phylogeny are most powerful. They also found that a significant source of error in divergence time estimation relates to misspecification of the clock model, especially when there are few informative calibrations. We cannot stress enough how sensitive posterior time estimates are to fossil information in constructing the time prior (Inoue et al. 2010). For example, different fossil calibration strategies have led to substantially different estimates of divergence times for mammals (dos Reis et al. 2014a) and early animals (dos Reis et al. 2015), regardless of how precise the branch length estimates are (Warnock et al. 2015). Thus, the field has reached the stage wherein there is general agreement that the choice of clock model and calibration strategy are fundamentally important to the accuracy of resulting age estimates, and thus, the way forward will clearly involve both empirical and simulation approaches to the problem. Here, we hope to contribute to this progress by conducting an exploration of model choice and calibration strategy in a classic empirical system: the Primates. Despite the fact that it is a relatively small and biologically uniform clade, primates have been inordinately and repeatedly the subject of divergence time analysis, with the first studies appearing at the very outset of molecular clock studies (Sarich and Wilson 1967), up to phylogenomic studies encompassing a large set of primate species (Perelman et al. 2011; Springer et al. 2012). This is largely due, undoubtedly, to the fact that we ourselves are members of this clade and can thus be forgiven for a persistent curiosity about our ancestral history. Age estimates for major primate divergence events have varied broadly among different studies (see Table 1), though one result has been relatively constant throughout: primate origins have been typically shown to predate the Cretaceous–Paleogene (K–Pg) mass-extinction event. Table 1. Overview of estimates of divergence times in Primates (in millions of years ago) for selected studies Study Data/Analysis Primates Haplorrhini Anthropoidea Platyrrhini Catarrhini Homininae Strepsirrhini Lorisiformes Lemuriformes Sarich and Wilson (1967) Immunological distance/strict clock 5 Hasegawa et al. (1985) 896 bp mtDNA/strict clock 3.7 Koop et al. (1989) B-globin DNA sequences/strict clock 55 40 25 Purvis (1995) Super Tree Analysis 57.2 39.9 14.7 8.1 41.8 22.1 39.6 Kumar and Hedges (1998) 658 nuclear genes (analyzed individually)/strict clock 47.6 6.7 Goodman et al. (1998) 60–80 kbp globin gene region/local clocks 63 58 40 7 50 23 45 Yang and Yoder (2003) 2404 bp mtDNA/local clocks 57.6 69.9 38.9 64.8 Poux and Douzery (2004) 1278 bp nDNA/local clocks 56.7, 58.4 45.4, 46.7 13.8, 14.2 39.6, 40.7 Eizirik et al. (2004) 8182 bp nDNA/relaxed clock 77.2 43.6 59.6 Steiper and Young (2006) 59.8 kbp of genomic data/ relaxed clock w/autocorrelated rates 77.5 42.9 20.8 30.5 8.6 57.1 Bininda-Emonds et al. (2007) 51,089 bp mtDNA and nDNA/ad hoc relaxed clock 87.7 Chatterjee et al. (2009) 6138 bp mtDNA and 2157 bp nDNA/relaxed clock 63.7 42.8 26.6 23.4 10.7 51.6 37.5 46.2 Perelman et al. (2011) 34,927 bp DNA/relaxed clock w/independent rates 87.2 81.3 43.5 24.8 31.6 8.3 68.7 40.3 58.6 Springer et al. (2012) 69 nDNA genes; 10 mtDNA genes/relaxed clock w/autocorrelated rates 67.8 61.2 40.6 23.3 25.1 8.0 54.2 34.7 50.0 dos Reis et al. (2012) 14,644 genes w/ 20.6 Mpb/relaxed clock w/autocorrelated rates 68.2 65.0 37.4 26.4 10.4 55.1 35.6 49.3 Finstermeier et al. (2013) complete mtDNA genomes/relaxed clock w/independent rates 66.2 63.0 45.3 22.0 32.0 8.4 56.9 34.5 47.1 Pozzi et al. (2014) complete mtDNA genomes/relaxed clock w/autocorrelated rates 74.1 70.0 46.7 20.9 32.1 10.6 66.3 40.3 59.6 This study (Strategy A) 3.4 million bp/relaxed clock w/ autocorrelated rates 74.4 70.6 45.0 25.3 32.6 10.5 62.7 37.9 57.2 This study (Strategy B) 3.4 million bp/relaxed clock w/ independent rates 84.8 78.6 45.9 25.3 28.8 7.6 64.0 38.2 55.3 Study Data/Analysis Primates Haplorrhini Anthropoidea Platyrrhini Catarrhini Homininae Strepsirrhini Lorisiformes Lemuriformes Sarich and Wilson (1967) Immunological distance/strict clock 5 Hasegawa et al. (1985) 896 bp mtDNA/strict clock 3.7 Koop et al. (1989) B-globin DNA sequences/strict clock 55 40 25 Purvis (1995) Super Tree Analysis 57.2 39.9 14.7 8.1 41.8 22.1 39.6 Kumar and Hedges (1998) 658 nuclear genes (analyzed individually)/strict clock 47.6 6.7 Goodman et al. (1998) 60–80 kbp globin gene region/local clocks 63 58 40 7 50 23 45 Yang and Yoder (2003) 2404 bp mtDNA/local clocks 57.6 69.9 38.9 64.8 Poux and Douzery (2004) 1278 bp nDNA/local clocks 56.7, 58.4 45.4, 46.7 13.8, 14.2 39.6, 40.7 Eizirik et al. (2004) 8182 bp nDNA/relaxed clock 77.2 43.6 59.6 Steiper and Young (2006) 59.8 kbp of genomic data/ relaxed clock w/autocorrelated rates 77.5 42.9 20.8 30.5 8.6 57.1 Bininda-Emonds et al. (2007) 51,089 bp mtDNA and nDNA/ad hoc relaxed clock 87.7 Chatterjee et al. (2009) 6138 bp mtDNA and 2157 bp nDNA/relaxed clock 63.7 42.8 26.6 23.4 10.7 51.6 37.5 46.2 Perelman et al. (2011) 34,927 bp DNA/relaxed clock w/independent rates 87.2 81.3 43.5 24.8 31.6 8.3 68.7 40.3 58.6 Springer et al. (2012) 69 nDNA genes; 10 mtDNA genes/relaxed clock w/autocorrelated rates 67.8 61.2 40.6 23.3 25.1 8.0 54.2 34.7 50.0 dos Reis et al. (2012) 14,644 genes w/ 20.6 Mpb/relaxed clock w/autocorrelated rates 68.2 65.0 37.4 26.4 10.4 55.1 35.6 49.3 Finstermeier et al. (2013) complete mtDNA genomes/relaxed clock w/independent rates 66.2 63.0 45.3 22.0 32.0 8.4 56.9 34.5 47.1 Pozzi et al. (2014) complete mtDNA genomes/relaxed clock w/autocorrelated rates 74.1 70.0 46.7 20.9 32.1 10.6 66.3 40.3 59.6 This study (Strategy A) 3.4 million bp/relaxed clock w/ autocorrelated rates 74.4 70.6 45.0 25.3 32.6 10.5 62.7 37.9 57.2 This study (Strategy B) 3.4 million bp/relaxed clock w/ independent rates 84.8 78.6 45.9 25.3 28.8 7.6 64.0 38.2 55.3 Notes: All estimates are mean estimates; see original works for confidence/credible intervals; all taxonomic designations signify crown nodes. Primates indicates the divergence of Strepsirrhini and Haplorrhini; Haplorrhini indicates the divergence of Tarsiidae from Anthropoidea; Platyrrhini; Catarrhini; Homininae indicates the divergence of Gorilla from Pan$$+$$Homo; Strepsirrhini indicates the divergence of Lorisiformes and Lemuriformes; Lorisiformes indicates the divergence of Galagidae and Lorisidae; Lemuriformes indicates the divergence of Daubentonia (the aye–aye) from other Malagasy lemurs. Table 1. Overview of estimates of divergence times in Primates (in millions of years ago) for selected studies Study Data/Analysis Primates Haplorrhini Anthropoidea Platyrrhini Catarrhini Homininae Strepsirrhini Lorisiformes Lemuriformes Sarich and Wilson (1967) Immunological distance/strict clock 5 Hasegawa et al. (1985) 896 bp mtDNA/strict clock 3.7 Koop et al. (1989) B-globin DNA sequences/strict clock 55 40 25 Purvis (1995) Super Tree Analysis 57.2 39.9 14.7 8.1 41.8 22.1 39.6 Kumar and Hedges (1998) 658 nuclear genes (analyzed individually)/strict clock 47.6 6.7 Goodman et al. (1998) 60–80 kbp globin gene region/local clocks 63 58 40 7 50 23 45 Yang and Yoder (2003) 2404 bp mtDNA/local clocks 57.6 69.9 38.9 64.8 Poux and Douzery (2004) 1278 bp nDNA/local clocks 56.7, 58.4 45.4, 46.7 13.8, 14.2 39.6, 40.7 Eizirik et al. (2004) 8182 bp nDNA/relaxed clock 77.2 43.6 59.6 Steiper and Young (2006) 59.8 kbp of genomic data/ relaxed clock w/autocorrelated rates 77.5 42.9 20.8 30.5 8.6 57.1 Bininda-Emonds et al. (2007) 51,089 bp mtDNA and nDNA/ad hoc relaxed clock 87.7 Chatterjee et al. (2009) 6138 bp mtDNA and 2157 bp nDNA/relaxed clock 63.7 42.8 26.6 23.4 10.7 51.6 37.5 46.2 Perelman et al. (2011) 34,927 bp DNA/relaxed clock w/independent rates 87.2 81.3 43.5 24.8 31.6 8.3 68.7 40.3 58.6 Springer et al. (2012) 69 nDNA genes; 10 mtDNA genes/relaxed clock w/autocorrelated rates 67.8 61.2 40.6 23.3 25.1 8.0 54.2 34.7 50.0 dos Reis et al. (2012) 14,644 genes w/ 20.6 Mpb/relaxed clock w/autocorrelated rates 68.2 65.0 37.4 26.4 10.4 55.1 35.6 49.3 Finstermeier et al. (2013) complete mtDNA genomes/relaxed clock w/independent rates 66.2 63.0 45.3 22.0 32.0 8.4 56.9 34.5 47.1 Pozzi et al. (2014) complete mtDNA genomes/relaxed clock w/autocorrelated rates 74.1 70.0 46.7 20.9 32.1 10.6 66.3 40.3 59.6 This study (Strategy A) 3.4 million bp/relaxed clock w/ autocorrelated rates 74.4 70.6 45.0 25.3 32.6 10.5 62.7 37.9 57.2 This study (Strategy B) 3.4 million bp/relaxed clock w/ independent rates 84.8 78.6 45.9 25.3 28.8 7.6 64.0 38.2 55.3 Study Data/Analysis Primates Haplorrhini Anthropoidea Platyrrhini Catarrhini Homininae Strepsirrhini Lorisiformes Lemuriformes Sarich and Wilson (1967) Immunological distance/strict clock 5 Hasegawa et al. (1985) 896 bp mtDNA/strict clock 3.7 Koop et al. (1989) B-globin DNA sequences/strict clock 55 40 25 Purvis (1995) Super Tree Analysis 57.2 39.9 14.7 8.1 41.8 22.1 39.6 Kumar and Hedges (1998) 658 nuclear genes (analyzed individually)/strict clock 47.6 6.7 Goodman et al. (1998) 60–80 kbp globin gene region/local clocks 63 58 40 7 50 23 45 Yang and Yoder (2003) 2404 bp mtDNA/local clocks 57.6 69.9 38.9 64.8 Poux and Douzery (2004) 1278 bp nDNA/local clocks 56.7, 58.4 45.4, 46.7 13.8, 14.2 39.6, 40.7 Eizirik et al. (2004) 8182 bp nDNA/relaxed clock 77.2 43.6 59.6 Steiper and Young (2006) 59.8 kbp of genomic data/ relaxed clock w/autocorrelated rates 77.5 42.9 20.8 30.5 8.6 57.1 Bininda-Emonds et al. (2007) 51,089 bp mtDNA and nDNA/ad hoc relaxed clock 87.7 Chatterjee et al. (2009) 6138 bp mtDNA and 2157 bp nDNA/relaxed clock 63.7 42.8 26.6 23.4 10.7 51.6 37.5 46.2 Perelman et al. (2011) 34,927 bp DNA/relaxed clock w/independent rates 87.2 81.3 43.5 24.8 31.6 8.3 68.7 40.3 58.6 Springer et al. (2012) 69 nDNA genes; 10 mtDNA genes/relaxed clock w/autocorrelated rates 67.8 61.2 40.6 23.3 25.1 8.0 54.2 34.7 50.0 dos Reis et al. (2012) 14,644 genes w/ 20.6 Mpb/relaxed clock w/autocorrelated rates 68.2 65.0 37.4 26.4 10.4 55.1 35.6 49.3 Finstermeier et al. (2013) complete mtDNA genomes/relaxed clock w/independent rates 66.2 63.0 45.3 22.0 32.0 8.4 56.9 34.5 47.1 Pozzi et al. (2014) complete mtDNA genomes/relaxed clock w/autocorrelated rates 74.1 70.0 46.7 20.9 32.1 10.6 66.3 40.3 59.6 This study (Strategy A) 3.4 million bp/relaxed clock w/ autocorrelated rates 74.4 70.6 45.0 25.3 32.6 10.5 62.7 37.9 57.2 This study (Strategy B) 3.4 million bp/relaxed clock w/ independent rates 84.8 78.6 45.9 25.3 28.8 7.6 64.0 38.2 55.3 Notes: All estimates are mean estimates; see original works for confidence/credible intervals; all taxonomic designations signify crown nodes. Primates indicates the divergence of Strepsirrhini and Haplorrhini; Haplorrhini indicates the divergence of Tarsiidae from Anthropoidea; Platyrrhini; Catarrhini; Homininae indicates the divergence of Gorilla from Pan$$+$$Homo; Strepsirrhini indicates the divergence of Lorisiformes and Lemuriformes; Lorisiformes indicates the divergence of Galagidae and Lorisidae; Lemuriformes indicates the divergence of Daubentonia (the aye–aye) from other Malagasy lemurs. Our study explores the effects of an autocorrelated versus an uncorrelated rate model on age estimates, and also explores the consequences of two different interpretations of both the age and the placement of key fossils with the living primate radiation. We apply these two strategies to a large phylogenomic data set for Primates (372 species and 3.4 million aligned base pairs). Until very recently, reliable calculations of branch lengths and age estimates within an analysis of this magnitude would have been beyond the capacity of computational methods. We have tackled many of these challenges by deploying the sequential Bayesian procedure used by dos Reis et al. (2012) wherein the posterior age estimates derived from a small taxonomic sample with genome-scale data are then deployed as priors for a subsequent analysis with many species and a much-reduced nucleotide sample. This procedure reduces the computational cost of a typical combined data analysis. It also helps to alleviate the concerns with the “missing data” problem as, in our approach, the sequence likelihood is only calculated for the species present in the alignment (Yang and Rannala 2006). The molecular timeline for primate evolution that emerges from this study can be interpreted with confidence. The data set is sufficiently large to provide highly precise branch length estimates, and the methods used are robust in accommodating violations of the molecular clock. The comparison of the two calibration strategies reveals their impact on the results by giving different age estimates within the tree, though the variation in inferred ages is not extreme. Which ages are considered most accurate will depend in large part on the degree of confidence in the fossils and their placement. As an unanticipated result of the study, the difference in age estimates for the deepest nodes of the phylogeny differ markedly when comparing the molecular rate models, with Bayesian model selection supporting the autocorrelated model. As with previous studies over the past several decades, the ancestral primate lineage is hypothesized to have survived the great K–Pg mass extinction event. Methods Bayesian estimates of divergence times of Primates were obtained using a supermatrix of molecular data with 372 species and 3.44 million base pairs (Mbp), combined with 17 fossil calibrations. The matrix is the result of merging the 372-species and 61 thousand base-pairs (kpb) data set of Springer et al. (2012) with a 10-species subset of the genome-scale alignment of dos Reis et al. (2012). Bayesian analyses were done using the program MCMCTree (Yang 2007). We assessed the robustness of time estimates by varying the clock model (strict clock, independent rates, and correlated rates) and by obtaining estimates under two fossil calibration strategies. Note that time estimates were obtained in two steps: in the first step, estimates were obtained for the small phylogeny of 10 species with a long alignment (3.38 Mbp). The marginal posterior of times was then used to construct the time prior in the second step for the 372-species phylogeny with a shorter alignment (61 kbp). This approach is almost the same as analysing the fully concatenated alignment in one step (3.38 Mbp $$+$$ 0.061 Mbp), but is computationally less demanding. All alignments, tree topology and fossil calibrations are available in Supplementary Material available on Dryad at http://dx.doi.org/10.5061/dryad.c020q. Sequence Alignment and Tree Topology Springer et al. (2012) alignment. We retrieved the sequence alignment of Springer et al. (2012), which is an extended version of the alignment of Perelman et al. (2011). The alignment has 372 species (367 primates and 5 outgroup species) and 79 gene segments (69 nuclear and 10 mitochondrial). The composite lagomorph sequence (an outgroup) was removed. We added a scandentian species (Tupaia belangeri), because its complete genome is available in the alignment of dos Reis et al. (2012), and because it has the 10 gene segments from the mitochondrial genome available (accession NC_002521). Many of the nuclear gene segments in the alignment of Springer et al. (2012) were mixtures of introns, exons, and UTRs, with out-of-frame indels in some exons. We manually curated the exons and separated the coding and noncoding segments of the alignment. These adjustments were necessary to facilitate an informed partition-based analysis of the data. Our modified version of Springer’s alignment was thus divided into six partitions: (i) 1st and 2nd codon positions for mitochondrial genes; (ii) 3rd positions for mitochondrial genes; (iii) mitochondrial RNA genes; (iv) 1st and 2nd codon positions for nuclear genes; (v) 3rd positions for nuclear genes; and (vi) noncoding segments of nuclear genes (UTRs and introns). The concatenated alignment has 372 species and is 61,132 base pairs long (Table 2). Our partitioning into codon positions and coding versus noncoding sequences follows established recommendations (Shapiro et al. 2006; Yang and Rannala 2006; Nascimento et al. 2017). Table 2. Sequence alignment summary Alignment Partition$$^{\mathrm{a}}$$ Sites Species Missing data$$^{\mathrm{b}}$$ Springer et al. 1. mit 1st$$+$$2nd 4816 330 61.4% 2. mit 3rd 2408 330 61.4% 3. mit RNA 2169 220 45.7% 4. nuclear 1st$$+$$2nd 16,309 239 53.8% 5. nuclear 3rd 8156 239 53.8% 6. nuclear noncoding 27,274 220 46.3% Partitions 1–6 61,132 372 51.1% (68.8%) dos Reis et al. 7. nuclear 1st$$+$$2nd 2,253,316 10 0.0% 8. nuclear 3rd 1,126,658 10 0.0% Partitions 7–8 3,379,974 10 0.0% (97.3%) Total 3,441,106 372 (96.8%) Alignment Partition$$^{\mathrm{a}}$$ Sites Species Missing data$$^{\mathrm{b}}$$ Springer et al. 1. mit 1st$$+$$2nd 4816 330 61.4% 2. mit 3rd 2408 330 61.4% 3. mit RNA 2169 220 45.7% 4. nuclear 1st$$+$$2nd 16,309 239 53.8% 5. nuclear 3rd 8156 239 53.8% 6. nuclear noncoding 27,274 220 46.3% Partitions 1–6 61,132 372 51.1% (68.8%) dos Reis et al. 7. nuclear 1st$$+$$2nd 2,253,316 10 0.0% 8. nuclear 3rd 1,126,658 10 0.0% Partitions 7–8 3,379,974 10 0.0% (97.3%) Total 3,441,106 372 (96.8%) $$^{\mathrm{a}}$$For topology estimation with RAxML, we used seven partitions: partitions 1 to 3, then 4 and 7 as one partition, 5 and 8 as one partition, and partition 6 divided into two: UTRs and Introns. $$^{\mathrm{b}}$$Numbers in brackets are the % missing data for the RAxML analysis. Note that MCMCTree only uses the species present in a partition to calculate the likelihood for the partition. In RAxML missing species in a partition are represented as sequences of only gaps in the partition, and thus the amount of missing data is larger. Table 2. Sequence alignment summary Alignment Partition$$^{\mathrm{a}}$$ Sites Species Missing data$$^{\mathrm{b}}$$ Springer et al. 1. mit 1st$$+$$2nd 4816 330 61.4% 2. mit 3rd 2408 330 61.4% 3. mit RNA 2169 220 45.7% 4. nuclear 1st$$+$$2nd 16,309 239 53.8% 5. nuclear 3rd 8156 239 53.8% 6. nuclear noncoding 27,274 220 46.3% Partitions 1–6 61,132 372 51.1% (68.8%) dos Reis et al. 7. nuclear 1st$$+$$2nd 2,253,316 10 0.0% 8. nuclear 3rd 1,126,658 10 0.0% Partitions 7–8 3,379,974 10 0.0% (97.3%) Total 3,441,106 372 (96.8%) Alignment Partition$$^{\mathrm{a}}$$ Sites Species Missing data$$^{\mathrm{b}}$$ Springer et al. 1. mit 1st$$+$$2nd 4816 330 61.4% 2. mit 3rd 2408 330 61.4% 3. mit RNA 2169 220 45.7% 4. nuclear 1st$$+$$2nd 16,309 239 53.8% 5. nuclear 3rd 8156 239 53.8% 6. nuclear noncoding 27,274 220 46.3% Partitions 1–6 61,132 372 51.1% (68.8%) dos Reis et al. 7. nuclear 1st$$+$$2nd 2,253,316 10 0.0% 8. nuclear 3rd 1,126,658 10 0.0% Partitions 7–8 3,379,974 10 0.0% (97.3%) Total 3,441,106 372 (96.8%) $$^{\mathrm{a}}$$For topology estimation with RAxML, we used seven partitions: partitions 1 to 3, then 4 and 7 as one partition, 5 and 8 as one partition, and partition 6 divided into two: UTRs and Introns. $$^{\mathrm{b}}$$Numbers in brackets are the % missing data for the RAxML analysis. Note that MCMCTree only uses the species present in a partition to calculate the likelihood for the partition. In RAxML missing species in a partition are represented as sequences of only gaps in the partition, and thus the amount of missing data is larger. dos Reis et al. (2012) alignment. We retrieved the genome-scale sequence alignment of dos Reis et al. (2012) of 36 mammal species, from which we extracted the sequences for 9 primates and 1 scandentian. The dos Reis et al. (2012) alignment was prepared using the highly curated mammalian genomes available in Ensembl. Though three additional primate genomes have become available in this database in the time since the original alignment was prepared, it is unlikely that their inclusion would change our results. The nine species represented in our study provide comprehensive phylogenetic representation of all major nodes in the primate tree and represent each of the higher-level clades (Fig. 1, inset). The original alignment has 14,632 nuclear, protein-coding genes, from which we removed 43 genes that were already present in the Springer alignment and 1 gene that was extremely long. All columns in the alignment with ambiguous nucleotides were removed, though care was taken not to disrupt the reading frame of the aligned coding sequences. The alignment was divided into two partitions: (i) 1st and 2nd codon positions; and (ii) 3rd codon positions. The final alignment has 10 species and is 3,441,106 base pairs long (missing data 0%, Table 2). Figure 1. View largeDownload slide The timetree of Primates. Nodes are drawn at their posterior mean ages in millions of years ago (Ma) estimated under the autocorrelated-rates (AR) clock model and calibration strategy A. Filled dots indicate nodes calibrated with the posterior times from the 10-species tree (inset figure), and empty dots indicate nodes with fossil constraints in the 372-species tree. Horizontal bars and numbers in parenthesis represent the 95% posterior CI for the node ages. Numbers associated with branches are ML bootstrap support values of major clades. Figure 1. View largeDownload slide The timetree of Primates. Nodes are drawn at their posterior mean ages in millions of years ago (Ma) estimated under the autocorrelated-rates (AR) clock model and calibration strategy A. Filled dots indicate nodes calibrated with the posterior times from the 10-species tree (inset figure), and empty dots indicate nodes with fossil constraints in the 372-species tree. Horizontal bars and numbers in parenthesis represent the 95% posterior CI for the node ages. Numbers associated with branches are ML bootstrap support values of major clades. Tree topology. The topology of the 372-species phylogeny was estimated by maximum likelihood (ML) using RAxML v 8.0.19 (Stamatakis 2014) under the GTR$$+$$G model (Yang 1994a,b), using seven partitions (Table 2) and 100 bootstrap replicates. Fossil Calibrations and Time Prior The two fossil calibration strategies used in this study are summarized in Table 3. They represent two different interpretations of the fossil record to construct calibrations for use in molecular clock dating analyses. Calibration strategy A is novel to this study, and calibration strategy B is based on the primate calibrations of dos Reis et al. (2012). Detailed justifications for the novel calibrations are provided in Appendix 1. Table 3. Fossil calibrations used in this study Calibration strategy$$^{\mathrm{a}}$$ Crown group Minimum (Ma) Maximum (Ma) MCMCTree Calibration$$^{\mathrm{b}}$$ Strategy A Human-Chimp$$^{\mathrm{c}}$$ 7.5 (†Sahelanthropus) 10 (unlikely before stem gorilla †Chororapithecus) B(0.075, 0.10, 0.01, 0.20) Human-Gorilla$$^{\mathrm{c}}$$ 10 (†Chororapithecus) 13.2 (unlikely before stem hominid †Sivapithecus) B(0.10, 0.132, 0.01, 0.20) Hominidae$$^{\mathrm{c}}$$ 11.2 (†Sivapithecus) 28 (unlikely before stem hominoid †Kamoyapithecus) B(0.112, 0.28, 0.01, 0.10) Catarrhini$$^{\mathrm{c}}$$ 25 (†Kamoyapithecus) 33.7 (absence of hominoids) B(0.25, 0.337, 0.01, 0.10) Anthropoidea$$^{\mathrm{c}}$$ 41 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 62.1 (K–Pg fossil modeling$$^{\mathrm{d}})$$ ST(0.4754, 0.0632, 0.98, 22.85) Strepsirrhini$$^{\mathrm{c}}$$ 37 (†Saharagalago) 58 (unlikely before †Altiatlasius) B(0.37, 0.58, 0.01, 0.10) Primates$$^{\mathrm{c}}$$ 57.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 88.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ S2N(0.474, 0.65, 0.0365, $$-$$3400, 0.6502, 0.1375, 11409) Euarchonta$$^{\mathrm{c}}$$ 65 (†Purgatorius) 130 (absence of placentals) G(36, 36.9) Strategy B Chimp-Human 5.7 (†Orrorin) 10 (absence of hominines) B(0.057, 0.10, 0.01, 0.05) Gorilla-Human 7.25 (†Chororapithecus) – L(0.0725, 0.1, 2) Hominidae 11.2 (†Sivapithecus) 33.7 (absence of pongines) B(0.112, 0.337, 0.05, 0.05) Catarrhini 23.5 (†Proconsul) 34 (absence of hominoids) B(0.235, 0.34, 0.01, 0.05) Anthropoidea 33.7 (†Catopithecus) — L(0.337, 0.1, 2) Haplorrrhini$$^{\mathrm{c}}$$ 45 (†Tarsius) — L(0.45, 0.1, 2) Strepsirrhini 33.7 (†Karanisia) 55.6 (absence of strepsirrhines) B(0.337, 0.556, 0.01, 0.05) Primates 55.6 (†Altiatlasius) — L(0.556, 0.1, 2) Euarchonta 61.5 (carpolestids and plesiadapids) 130 (absence of placentals) B(0.615, 1.30, 0.01, 0.05) Lorises$$^{\mathrm{c}}$$ 14 (Lorisidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.14, 0.37, 0.01, 0.10) Galagos$$^{\mathrm{c}}$$ 15 (Galagidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.15, 0.37, 0.01, 0.10) Lorisiformes$$^{\mathrm{c}}$$ 18 (†Mioeouticus) 38 (unlikely before †Karanisia/†Saharagalago) B(0.18, 0.38, 0.01, 0.10) Shared by Platyrrhini$$^{\mathrm{c}}$$ 15.7 (stem Pitheciinae) 33 (unlikely before †Catopithecus) B(0.157, 0.33, 0.01, 0.10) both Atelidae$$^{\mathrm{c}}$$ 12.8 (†Stirtonia) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) strategies Cebidae$$^{\mathrm{c}}$$ 12.8 (†Neosaimiri) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) Cercopithecinae$$^{\mathrm{c}}$$ 5 (†Parapapio) 23 (unlikely before †Prohylobates/Kamoyapithecus) B(0.05, 0.23, 0.01, 0.10) Colobinae$$^{\mathrm{c}}$$ 9.8 (†Microcolobus) 23 (unlikely before †Prohylobates/Kamoyapotheics) B(0.098, 0.23, 0.01, 0.10) Calibration strategy$$^{\mathrm{a}}$$ Crown group Minimum (Ma) Maximum (Ma) MCMCTree Calibration$$^{\mathrm{b}}$$ Strategy A Human-Chimp$$^{\mathrm{c}}$$ 7.5 (†Sahelanthropus) 10 (unlikely before stem gorilla †Chororapithecus) B(0.075, 0.10, 0.01, 0.20) Human-Gorilla$$^{\mathrm{c}}$$ 10 (†Chororapithecus) 13.2 (unlikely before stem hominid †Sivapithecus) B(0.10, 0.132, 0.01, 0.20) Hominidae$$^{\mathrm{c}}$$ 11.2 (†Sivapithecus) 28 (unlikely before stem hominoid †Kamoyapithecus) B(0.112, 0.28, 0.01, 0.10) Catarrhini$$^{\mathrm{c}}$$ 25 (†Kamoyapithecus) 33.7 (absence of hominoids) B(0.25, 0.337, 0.01, 0.10) Anthropoidea$$^{\mathrm{c}}$$ 41 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 62.1 (K–Pg fossil modeling$$^{\mathrm{d}})$$ ST(0.4754, 0.0632, 0.98, 22.85) Strepsirrhini$$^{\mathrm{c}}$$ 37 (†Saharagalago) 58 (unlikely before †Altiatlasius) B(0.37, 0.58, 0.01, 0.10) Primates$$^{\mathrm{c}}$$ 57.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 88.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ S2N(0.474, 0.65, 0.0365, $$-$$3400, 0.6502, 0.1375, 11409) Euarchonta$$^{\mathrm{c}}$$ 65 (†Purgatorius) 130 (absence of placentals) G(36, 36.9) Strategy B Chimp-Human 5.7 (†Orrorin) 10 (absence of hominines) B(0.057, 0.10, 0.01, 0.05) Gorilla-Human 7.25 (†Chororapithecus) – L(0.0725, 0.1, 2) Hominidae 11.2 (†Sivapithecus) 33.7 (absence of pongines) B(0.112, 0.337, 0.05, 0.05) Catarrhini 23.5 (†Proconsul) 34 (absence of hominoids) B(0.235, 0.34, 0.01, 0.05) Anthropoidea 33.7 (†Catopithecus) — L(0.337, 0.1, 2) Haplorrrhini$$^{\mathrm{c}}$$ 45 (†Tarsius) — L(0.45, 0.1, 2) Strepsirrhini 33.7 (†Karanisia) 55.6 (absence of strepsirrhines) B(0.337, 0.556, 0.01, 0.05) Primates 55.6 (†Altiatlasius) — L(0.556, 0.1, 2) Euarchonta 61.5 (carpolestids and plesiadapids) 130 (absence of placentals) B(0.615, 1.30, 0.01, 0.05) Lorises$$^{\mathrm{c}}$$ 14 (Lorisidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.14, 0.37, 0.01, 0.10) Galagos$$^{\mathrm{c}}$$ 15 (Galagidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.15, 0.37, 0.01, 0.10) Lorisiformes$$^{\mathrm{c}}$$ 18 (†Mioeouticus) 38 (unlikely before †Karanisia/†Saharagalago) B(0.18, 0.38, 0.01, 0.10) Shared by Platyrrhini$$^{\mathrm{c}}$$ 15.7 (stem Pitheciinae) 33 (unlikely before †Catopithecus) B(0.157, 0.33, 0.01, 0.10) both Atelidae$$^{\mathrm{c}}$$ 12.8 (†Stirtonia) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) strategies Cebidae$$^{\mathrm{c}}$$ 12.8 (†Neosaimiri) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) Cercopithecinae$$^{\mathrm{c}}$$ 5 (†Parapapio) 23 (unlikely before †Prohylobates/Kamoyapithecus) B(0.05, 0.23, 0.01, 0.10) Colobinae$$^{\mathrm{c}}$$ 9.8 (†Microcolobus) 23 (unlikely before †Prohylobates/Kamoyapotheics) B(0.098, 0.23, 0.01, 0.10) $$^{\mathrm{a}}$$Calibration strategies A and B are applied to the nodes in the phylogeny of 10 species. The shared calibrations are applied to the large tree of 372 species. Ages are in millions of years ago (Ma). $$^{\mathrm{b}}$$B($$t_{\mathrm{L}}, t_{\mathrm{U}}, p_{\mathrm{L}}, p_{\mathrm{U}})$$ means the node age is calibrated by a uniform distribution bounded between a minimum time $$t_{\mathrm{L}}$$, and a maximum time $$t_{\mathrm{U}}$$, with probabilities $$p_{\mathrm{L}}$$ and $$p_{\mathrm{U}}$$ that the age is outside the bounds. L($$t_{\mathrm{L}}$$, $$p$$, $$c$$, $$p_{\mathrm{L}})$$ means the node age is calibrated by a truncated Cauchy distribution with minimum age $$t_{\mathrm{L}}$$ and parameters $$p$$ and $$c$$, with the probability that the age is younger than the minimum bound to be $$p_{\mathrm{L}} \quad =$$ 5% (Inoue et al. 2010). ST($$a$$, $$b$$, $$c$$, $$d)$$ means the node age is calibrated by a skew-t density with parameters $$a$$, $$b$$, $$c,$$ and $$d$$ (Wilkinson et al. 2011). S2N($$a$$, $$b$$, $$c)$$ means the node age is calibrated by a mixture of two skew-normal distributions (Wilkinson et al. 2011). G($$a$$, $$b)$$ means the node age is calibrated by a gamma distribution with shape $$a$$ and rate $$b$$. $$^{\mathrm{c}}$$Detailed justifications for these fossil calibrations are given in the Appendix. For all other calibrations, justifications are in Benton et al. (2009) and dos Reis et al. (2012). $$^{\mathrm{d}}$$Calibration densities are the posterior distribution from a model of fossil preservation and discovery with species diversification that takes into account the effects of the K–Pg extinction in the model (Wilkinson et al. 2011). Table 3. Fossil calibrations used in this study Calibration strategy$$^{\mathrm{a}}$$ Crown group Minimum (Ma) Maximum (Ma) MCMCTree Calibration$$^{\mathrm{b}}$$ Strategy A Human-Chimp$$^{\mathrm{c}}$$ 7.5 (†Sahelanthropus) 10 (unlikely before stem gorilla †Chororapithecus) B(0.075, 0.10, 0.01, 0.20) Human-Gorilla$$^{\mathrm{c}}$$ 10 (†Chororapithecus) 13.2 (unlikely before stem hominid †Sivapithecus) B(0.10, 0.132, 0.01, 0.20) Hominidae$$^{\mathrm{c}}$$ 11.2 (†Sivapithecus) 28 (unlikely before stem hominoid †Kamoyapithecus) B(0.112, 0.28, 0.01, 0.10) Catarrhini$$^{\mathrm{c}}$$ 25 (†Kamoyapithecus) 33.7 (absence of hominoids) B(0.25, 0.337, 0.01, 0.10) Anthropoidea$$^{\mathrm{c}}$$ 41 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 62.1 (K–Pg fossil modeling$$^{\mathrm{d}})$$ ST(0.4754, 0.0632, 0.98, 22.85) Strepsirrhini$$^{\mathrm{c}}$$ 37 (†Saharagalago) 58 (unlikely before †Altiatlasius) B(0.37, 0.58, 0.01, 0.10) Primates$$^{\mathrm{c}}$$ 57.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 88.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ S2N(0.474, 0.65, 0.0365, $$-$$3400, 0.6502, 0.1375, 11409) Euarchonta$$^{\mathrm{c}}$$ 65 (†Purgatorius) 130 (absence of placentals) G(36, 36.9) Strategy B Chimp-Human 5.7 (†Orrorin) 10 (absence of hominines) B(0.057, 0.10, 0.01, 0.05) Gorilla-Human 7.25 (†Chororapithecus) – L(0.0725, 0.1, 2) Hominidae 11.2 (†Sivapithecus) 33.7 (absence of pongines) B(0.112, 0.337, 0.05, 0.05) Catarrhini 23.5 (†Proconsul) 34 (absence of hominoids) B(0.235, 0.34, 0.01, 0.05) Anthropoidea 33.7 (†Catopithecus) — L(0.337, 0.1, 2) Haplorrrhini$$^{\mathrm{c}}$$ 45 (†Tarsius) — L(0.45, 0.1, 2) Strepsirrhini 33.7 (†Karanisia) 55.6 (absence of strepsirrhines) B(0.337, 0.556, 0.01, 0.05) Primates 55.6 (†Altiatlasius) — L(0.556, 0.1, 2) Euarchonta 61.5 (carpolestids and plesiadapids) 130 (absence of placentals) B(0.615, 1.30, 0.01, 0.05) Lorises$$^{\mathrm{c}}$$ 14 (Lorisidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.14, 0.37, 0.01, 0.10) Galagos$$^{\mathrm{c}}$$ 15 (Galagidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.15, 0.37, 0.01, 0.10) Lorisiformes$$^{\mathrm{c}}$$ 18 (†Mioeouticus) 38 (unlikely before †Karanisia/†Saharagalago) B(0.18, 0.38, 0.01, 0.10) Shared by Platyrrhini$$^{\mathrm{c}}$$ 15.7 (stem Pitheciinae) 33 (unlikely before †Catopithecus) B(0.157, 0.33, 0.01, 0.10) both Atelidae$$^{\mathrm{c}}$$ 12.8 (†Stirtonia) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) strategies Cebidae$$^{\mathrm{c}}$$ 12.8 (†Neosaimiri) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) Cercopithecinae$$^{\mathrm{c}}$$ 5 (†Parapapio) 23 (unlikely before †Prohylobates/Kamoyapithecus) B(0.05, 0.23, 0.01, 0.10) Colobinae$$^{\mathrm{c}}$$ 9.8 (†Microcolobus) 23 (unlikely before †Prohylobates/Kamoyapotheics) B(0.098, 0.23, 0.01, 0.10) Calibration strategy$$^{\mathrm{a}}$$ Crown group Minimum (Ma) Maximum (Ma) MCMCTree Calibration$$^{\mathrm{b}}$$ Strategy A Human-Chimp$$^{\mathrm{c}}$$ 7.5 (†Sahelanthropus) 10 (unlikely before stem gorilla †Chororapithecus) B(0.075, 0.10, 0.01, 0.20) Human-Gorilla$$^{\mathrm{c}}$$ 10 (†Chororapithecus) 13.2 (unlikely before stem hominid †Sivapithecus) B(0.10, 0.132, 0.01, 0.20) Hominidae$$^{\mathrm{c}}$$ 11.2 (†Sivapithecus) 28 (unlikely before stem hominoid †Kamoyapithecus) B(0.112, 0.28, 0.01, 0.10) Catarrhini$$^{\mathrm{c}}$$ 25 (†Kamoyapithecus) 33.7 (absence of hominoids) B(0.25, 0.337, 0.01, 0.10) Anthropoidea$$^{\mathrm{c}}$$ 41 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 62.1 (K–Pg fossil modeling$$^{\mathrm{d}})$$ ST(0.4754, 0.0632, 0.98, 22.85) Strepsirrhini$$^{\mathrm{c}}$$ 37 (†Saharagalago) 58 (unlikely before †Altiatlasius) B(0.37, 0.58, 0.01, 0.10) Primates$$^{\mathrm{c}}$$ 57.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ 88.6 (K–Pg fossil modeling$$^{\mathrm{d}})$$ S2N(0.474, 0.65, 0.0365, $$-$$3400, 0.6502, 0.1375, 11409) Euarchonta$$^{\mathrm{c}}$$ 65 (†Purgatorius) 130 (absence of placentals) G(36, 36.9) Strategy B Chimp-Human 5.7 (†Orrorin) 10 (absence of hominines) B(0.057, 0.10, 0.01, 0.05) Gorilla-Human 7.25 (†Chororapithecus) – L(0.0725, 0.1, 2) Hominidae 11.2 (†Sivapithecus) 33.7 (absence of pongines) B(0.112, 0.337, 0.05, 0.05) Catarrhini 23.5 (†Proconsul) 34 (absence of hominoids) B(0.235, 0.34, 0.01, 0.05) Anthropoidea 33.7 (†Catopithecus) — L(0.337, 0.1, 2) Haplorrrhini$$^{\mathrm{c}}$$ 45 (†Tarsius) — L(0.45, 0.1, 2) Strepsirrhini 33.7 (†Karanisia) 55.6 (absence of strepsirrhines) B(0.337, 0.556, 0.01, 0.05) Primates 55.6 (†Altiatlasius) — L(0.556, 0.1, 2) Euarchonta 61.5 (carpolestids and plesiadapids) 130 (absence of placentals) B(0.615, 1.30, 0.01, 0.05) Lorises$$^{\mathrm{c}}$$ 14 (Lorisidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.14, 0.37, 0.01, 0.10) Galagos$$^{\mathrm{c}}$$ 15 (Galagidae gen et sp. nov) 37 (unlikely before †Karanisia) B(0.15, 0.37, 0.01, 0.10) Lorisiformes$$^{\mathrm{c}}$$ 18 (†Mioeouticus) 38 (unlikely before †Karanisia/†Saharagalago) B(0.18, 0.38, 0.01, 0.10) Shared by Platyrrhini$$^{\mathrm{c}}$$ 15.7 (stem Pitheciinae) 33 (unlikely before †Catopithecus) B(0.157, 0.33, 0.01, 0.10) both Atelidae$$^{\mathrm{c}}$$ 12.8 (†Stirtonia) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) strategies Cebidae$$^{\mathrm{c}}$$ 12.8 (†Neosaimiri) 18 (unlikely before †Soriacebus) B(0.128, 0.18, 0.01, 0.10) Cercopithecinae$$^{\mathrm{c}}$$ 5 (†Parapapio) 23 (unlikely before †Prohylobates/Kamoyapithecus) B(0.05, 0.23, 0.01, 0.10) Colobinae$$^{\mathrm{c}}$$ 9.8 (†Microcolobus) 23 (unlikely before †Prohylobates/Kamoyapotheics) B(0.098, 0.23, 0.01, 0.10) $$^{\mathrm{a}}$$Calibration strategies A and B are applied to the nodes in the phylogeny of 10 species. The shared calibrations are applied to the large tree of 372 species. Ages are in millions of years ago (Ma). $$^{\mathrm{b}}$$B($$t_{\mathrm{L}}, t_{\mathrm{U}}, p_{\mathrm{L}}, p_{\mathrm{U}})$$ means the node age is calibrated by a uniform distribution bounded between a minimum time $$t_{\mathrm{L}}$$, and a maximum time $$t_{\mathrm{U}}$$, with probabilities $$p_{\mathrm{L}}$$ and $$p_{\mathrm{U}}$$ that the age is outside the bounds. L($$t_{\mathrm{L}}$$, $$p$$, $$c$$, $$p_{\mathrm{L}})$$ means the node age is calibrated by a truncated Cauchy distribution with minimum age $$t_{\mathrm{L}}$$ and parameters $$p$$ and $$c$$, with the probability that the age is younger than the minimum bound to be $$p_{\mathrm{L}} \quad =$$ 5% (Inoue et al. 2010). ST($$a$$, $$b$$, $$c$$, $$d)$$ means the node age is calibrated by a skew-t density with parameters $$a$$, $$b$$, $$c,$$ and $$d$$ (Wilkinson et al. 2011). S2N($$a$$, $$b$$, $$c)$$ means the node age is calibrated by a mixture of two skew-normal distributions (Wilkinson et al. 2011). G($$a$$, $$b)$$ means the node age is calibrated by a gamma distribution with shape $$a$$ and rate $$b$$. $$^{\mathrm{c}}$$Detailed justifications for these fossil calibrations are given in the Appendix. For all other calibrations, justifications are in Benton et al. (2009) and dos Reis et al. (2012). $$^{\mathrm{d}}$$Calibration densities are the posterior distribution from a model of fossil preservation and discovery with species diversification that takes into account the effects of the K–Pg extinction in the model (Wilkinson et al. 2011). Fossil calibration strategy A. We used the fossil-based prior densities constructed by Wilkinson et al. (2011) to calibrate the ages of crown Primates and crown Anthropoidea. The prior densities were constructed by modeling the processes of speciation, extinction, fossil preservation, and rates of fossil discovery in Primates. The effects of the K–Pg extinction were accounted for in the model. We calibrated six more node ages by using uniform distribution densities with soft bounds (Yang and Rannala 2006). We set the probability of violating a minimum bound to 1%. Because maximum bounds are based on weak evidence, we set the probability that a maximum bound is violated to 10% or 20%. The crown Haplorrhini node was left with no calibration as the branch separating that clade from crown Primates is very short and we wanted to avoid truncation with the fossil-modeling density on crown Primates. The prior on the age of crown Haplorrhini is instead set using the birth–death process with parameters $$\lambda =\mu =1$$ and $$\rho =0$$. These parameter values specify a uniform kernel density (Yang and Rannala 1997, equation 7). Fossil calibration strategy B. We used the same nine calibrations that dos Reis et al. (2012) used to calibrate the Primates and Scandentia clades. An additional calibration based on †Tarsius sp$$.$$ (Beard et al. 1994) was used for the Haplorrhini node. For nodes with a minimum bound only, modeled using a truncated Cauchy density, the spread parameter was set to $$c=2$$ (Inoue et al. 2010). For maximum bounds the probability that the bound was violated was set to 5%. There are other differences between strategies A and B (Table 1). For example, in A, we considered †Sahelanthropus, dated to 7.25 million years ago (Ma), to be the oldest member of the human-chimpanzee clade (Brunet et al. 2002, 2005) and used it to calibrate the clade accordingly, while in B, dos Reis et al. (2012) used †Orrorin (5.7 Ma) instead. In A, †Chororapithecus is given an age of 10 Ma (Geraads et al. 2002), while in B it is given the younger (perhaps more conservative) age of 7.25 Ma (Benton et al. 2009). We note that the ages of fossils and their relationships to extant groups are often controversial, and we cannot overemphasize the degree to which differences of opinion among palaeontologists are an important source of uncertainty in the construction of fossil calibrations, and accordingly, divergence time estimates throughout the phylogeny. Calibrating the 372-species phylogeny. Strategies A and B were used to obtain time estimates for the 10-species phylogeny using the 3.38 Mbp alignment. Then skew-t densities were fitted by ML to the marginal posterior ages of each of the 9 internal nodes in the 10-species phylogeny and used to calibrate the corresponding nodes in the 372-species tree. Eight additional fossil calibrations (Table 3) were used to calibrate additional nodes in the 372-species tree. For nodes without calibrations, the time prior was constructed using the birth–death process with parameters $$\lambda =\mu =1$$ and $$\rho =0$$. Bayesian time estimation then proceeded on the 372-species tree and 61 kbp alignment as usual. Rate Prior For the 10-species analysis, the rate prior was set as follows: the nuclear substitution rate at third codon positions in apes is roughly within 10$$^{\mathrm{-9}}$$ substitutions per site per year (s/s/y) (Burgess and Yang 2008). At first and second codon positions, it is about a fifth of the third position rate, or $$2\times 10^{\mathrm{-10}}$$ s/s/y. This gives roughly an overall rate of about $$5\times 10^{\mathrm{-10}}$$ s/s/y for the three positions combined. We thus used a diffuse gamma density G(2, 40) with mean 0.05 and 95% prior credibility interval (CI) 0.00606–0.139 (our time unit is 100 My, thus, this corresponds to $$6.06\times 10^{\mathrm{-12}}$$ to $$1.39\times 10^{\mathrm{-10}}$$ s/s/y). The analysis was conducted under both the autocorrelated rates (AR) and independent rates (IR) models. Parameter $$\sigma^{\mathrm{2}}$$ in the AR and IR models was assigned a gamma prior G(1, 10). Note that the average rate for loci, $$\mu_{i}$$, and $$\sigma_{i}^{2}$$ are assigned a gamma-Dirichlet prior (dos Reis et al. 2014b). For the 372-species phylogeny, the rate prior was assigned as follows: the mitochondrial substitution rate at third positions is about 20 times the rate at third positions in nuclear genes or $$2\times 10^{\mathrm{-8}}$$. Assuming 1st and 2nd codon positions evolve at about a fifth of the third position rate we get roughly $$4\times 10^{\mathrm{-9}}$$. The prior mean is then approximately $$2.5\times 10^{\mathrm{-9}}$$ s/s/y, which is the weighted average (by number of sites) of the substitution rates for the nuclear and mitochondrial partitions. We thus used a gamma density G(2, 8) with mean 0.25 and 95% CI 0.0302–0.696. For $$\sigma^{\mathrm{2}}$$ we used G(1, 10). MCMC and Bayesian Selection of Clock Model MCMC analyses were carried out with the program MCMCTree (Yang 2007), using the approximate likelihood method (dos Reis and Yang 2011). Convergence of the MCMC to the posterior distribution was assessed by running the analyses multiple times. MCMCtree runs were carried out without sequence data to calculate the joint prior of node ages. Results from all analyses were summarized as posterior means and 95% CIs. We have implemented marginal likelihood calculation by thermodynamic integration (path sampling) in the program MCMCTree. This allows us to calculate Bayes factors (BF) and posterior model probabilities to select for a clock model in the analysis. Details of our implementation are given in Appendix 2. Extensive discussions on marginal likelihood estimation by thermodynamic integration and stepping-stones (a related method) in the phylogenetic context are given in Lartillot and Philippe (2006), Lepage et al. (2007), Xie et al. (2011), and Baele et al. (2012a,b). A detailed simulation study is given in Ho et al. (2015). Thermodynamic integration is computationally intensive as we must sample from the power posterior $$f(\theta) f(D$$| $$\theta )^{\beta }$$, in a sampling path from the prior ($$\beta =0$$) to the posterior ($$\beta=1$$). Because the approximation to the likelihood is not good when samples are taken far away from the ML estimate (dos Reis and Yang 2011), as it happens when $$\beta$$ is small, the approximation cannot be used in the calculation of the power posterior. Thus, we use exact likelihood calculation on a smaller data set of nine primate species (inset of Fig. 1), for the six partitions of the Springer alignment (Table 2) to perform the Bayesian selection of clock model. We use 64 $$\quad$$$$\beta$$-points to construct the sampling path from the prior to the posterior and calculate the marginal likelihoods for the strict clock, and the AR and IR models. Effect of Genome-Scale Data In a conventional statistical inference problem, the variance of an estimate decreases in proportion to 1/$$n$$, with $$n$$ to be the sample size. Thus, as the sample size approaches infinity, the variance of an estimate approaches zero and the estimate converges to the true value. In divergence time estimation, which is an unconventional estimation problem, the nonidentifiability of times and rates means that the uncertainty in the posterior of times does not converge to zero as the amount of molecular data (the sample size) approaches infinity, but rather converge to a limiting value imposed by the uncertainties in the fossil calibrations (Yang and Rannala 2006; Rannala and Yang 2007). For infinitely long alignments, an infinite-sites plot (a plot of uncertainty in the time posterior, measured as the width of the CI, i.e., the difference between the 2.5% and 97.5% limits vs. the mean posterior of time) would converge onto a straight line. This line represents the amount of uncertainty in time estimates for every 1 million years (My) of divergence that is due solely to uncertainties in the fossil calibrations. We calculate the infinite-sites plots for time estimates on the 372-species phylogeny to study the effect of genome-scale data on the uncertainty of species-level time estimates. Results A Timeline of Primate Evolution Figure 1 illustrates time estimates under the AR model and calibration strategy A (Figs. 2–4 show detailed timetrees for the major clades). Under calibration strategy A and the AR model, we find that crown Primates originated 79.2–70.0 million years ago (Ma), before the K–Pg event at 66 Ma. However, the diversification of the main clades occurred much later. Crown Anthropoidea originated 48.3–41.8 Ma, with its two main crown groups, Catarrhini (Old World monkeys and apes) and Platyrrhini (New World monkeys) originating at 35.1–30.4 Ma and 27.5–23.6 Ma, respectively. Crown tarsiers originated 33.5–15.5 Ma. Crown Strepsirrhini date back to 66.8–58.8 Ma, with its two main crown groups, Lemuriformes and Lorisiformes, dating back to 61.6–52.7 Ma and 40.9–34.1 Ma, respectively. Figure 2. View largeDownload slide Strepsirrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 2. View largeDownload slide Strepsirrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 3. View largeDownload slide Catarrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 3. View largeDownload slide Catarrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 4. View largeDownload slide Tarsiidae and platyrrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Figure 4. View largeDownload slide Tarsiidae and platyrrhine portion of the primate timetree (AR clock and calibration strategy A). Legend as for Figure 1. Calibration strategy B under the AR model gives similar node age estimates for the younger nodes in the tree (i.e., the 95% CI of node age overlap, Fig. 5A). However, for the older nodes in the phylogeny (and in particular for Euarchonta, Primatomorpha, Primates, Haplorrhini, Lemuriformes, and Lemuriformes minus aye-aye), strategy A produced older estimates (Fig. 5A). Under strategy B a pre-K–Pg origin of crown Primates is also favored, although the posterior distribution of the age of crown Primates straddles the K–Pg boundary (71.4–63.9 Ma). The posterior probability for a pre-K–Pg origin of crown Primates is 80.0% under strategy B and 100% under strategy A. Posterior time estimates for all nodes under both strategies are given in the Supplementary Material spreadsheet available on Dryad. Figure 5. View largeDownload slide Effect of calibration strategy and relaxed-clock model. A) Posterior time estimates under fossil calibration strategy A versus time estimates under strategy B, for the AR clock model. B) Posterior time estimates under the AR versus IR clock models, for calibration strategy A. Figure 5. View largeDownload slide Effect of calibration strategy and relaxed-clock model. A) Posterior time estimates under fossil calibration strategy A versus time estimates under strategy B, for the AR clock model. B) Posterior time estimates under the AR versus IR clock models, for calibration strategy A. Note that the two calibration strategies are in many cases based on the same fossils (Table 3), and the intervals defined by the fossil bounds overlap extensively between the two strategies. However, the seemingly small differences between the two strategies lead to noticeable differences in the posterior time estimates (Fig. 5A). In general, minimum bound constraints are older in strategy A than in strategy B (Table 3), and thus this may be the cause of the older time estimates in A versus B. Time prior and effect of truncation and outgroups. User-specified calibration densities usually do not satisfy the constraint that descendant nodes must be younger than their ancestors, thus the dating methodology must “truncate” the calibration densities to satisfy the constraint to construct the time prior (Rannala 2016). The result is that user-specified calibration densities and marginal priors of node ages may look substantially different. Supplementary Figure S1 available on Dryad illustrates the effect of truncation on prior densities for strategies A and B. For example, in strategy B, the calibration densities on Euarchonta (the root of the phylogeny) and on Primates interact (the primate node has a Cauchy calibration density with a heavy tail), and consequently the prior density on the age of Euarchonta is pushed back (Supplementary Figure S1 available on Dryad). The result is that the marginal prior age of Euarchonta ranges from 136–178 Ma (Supplementary Figure S1 available on Dryad) instead of 130–61.5 Ma as in the calibration density (Table 3), while the upper age for the Primate prior is too old (127 Ma). In contrast, under strategy A, the calibration density on Primates has a much lighter tail, and thus the truncation effect with the Euarchonta node is minimal. The result is that the marginal time prior and the corresponding calibration densities for the Primates and Euarchonta nodes are very similar (Supplementary Figure S1 available on Dryad). Similarly, under strategy B, the priors for two other nodes (Anthropoidea and Human-Gorilla) that use the heavy-tailed Cauchy calibrations have upper 95% limits that also appear unreasonably old (86.3 Ma and 25.0 Ma, respectively). In general, calibration strategy A, which avoids using the long-tailed Cauchy calibrations, has calibration densities that are much closer to the resulting marginal priors, and thus strategy A results in a time prior which is much closer to the fossil information as interpreted by the palaeontologist. A set of calibrations for mitogenomic phylogenetic analysis. Mitogenomic markers are widely used to construct phylogenies of closely related primate species with examples seen in phylogeographic studies of diversification of primates in the Amazon (Nascimento et al. 2014) and in the timing of human diversification (Rieux et al. 2014). The posterior distributions obtained here for the 10-species genomic data are useful calibrations for mitogenomic studies. Note that these cannot be used if the molecular alignment contains nuclear data as the calibrations already contain the information from nuclear genomes. The list of skew-t calibrations is provided in Table 4, together with approximations based on the gamma distribution, which can be used in software that does not implement skew-t calibrations (such as BEAST or MrBayes). Table 4. Suggested skew-t and gamma calibrations for mitogenomic studies Crown group Skew-t$$^{\mathrm{a,b}}$$ Gamma$$^{\mathrm{a,c}}$$ Primates ST(0.878, 0.169, 2.41, 94.6) G(78.6, 77.6) Strepsirrhini ST(0.580, 0.0567, 2.12, 149) G(262, 327) Haplorrhini ST(0.705, 0.0628, 1.66, 960) G(254, 411) Anthropoidea ST(0.415, 0.0291, 0.949, 294) G(271, 363) Catarrhini ST(0.292, 0.0206, 0.995, 167) G(316, 733) Hominoidea ST(0.185, 0.0167, 2.44, 312) G(316, 1040) Human-Gorilla ST(0.0996, 0.0103, 19.6, 100) G(311, 1575) Human-Chimp ST(0.0788, 0.00687, 3.51, 6.15) G(292, 2715) Crown group Skew-t$$^{\mathrm{a,b}}$$ Gamma$$^{\mathrm{a,c}}$$ Primates ST(0.878, 0.169, 2.41, 94.6) G(78.6, 77.6) Strepsirrhini ST(0.580, 0.0567, 2.12, 149) G(262, 327) Haplorrhini ST(0.705, 0.0628, 1.66, 960) G(254, 411) Anthropoidea ST(0.415, 0.0291, 0.949, 294) G(271, 363) Catarrhini ST(0.292, 0.0206, 0.995, 167) G(316, 733) Hominoidea ST(0.185, 0.0167, 2.44, 312) G(316, 1040) Human-Gorilla ST(0.0996, 0.0103, 19.6, 100) G(311, 1575) Human-Chimp ST(0.0788, 0.00687, 3.51, 6.15) G(292, 2715) $$^{\mathrm{a}}$$Densities calculated under fossil calibration strategy A and AR model using a time unit of 100 My. $$^{\mathrm{b}}$$The parameters of the skew-t distribution are location, scale, shape, and df. $$^{\mathrm{c}}$$Note that here we use the shape ($$\alpha )$$ and rate ($$\beta )$$ parameterization. For the scale parameterization use $$s=1$$/$$\beta$$. The mean is $$\alpha$$/$$\beta$$ and variance is $$\alpha$$/$$\beta^{\mathrm{2}}$$. Table 4. Suggested skew-t and gamma calibrations for mitogenomic studies Crown group Skew-t$$^{\mathrm{a,b}}$$ Gamma$$^{\mathrm{a,c}}$$ Primates ST(0.878, 0.169, 2.41, 94.6) G(78.6, 77.6) Strepsirrhini ST(0.580, 0.0567, 2.12, 149) G(262, 327) Haplorrhini ST(0.705, 0.0628, 1.66, 960) G(254, 411) Anthropoidea ST(0.415, 0.0291, 0.949, 294) G(271, 363) Catarrhini ST(0.292, 0.0206, 0.995, 167) G(316, 733) Hominoidea ST(0.185, 0.0167, 2.44, 312) G(316, 1040) Human-Gorilla ST(0.0996, 0.0103, 19.6, 100) G(311, 1575) Human-Chimp ST(0.0788, 0.00687, 3.51, 6.15) G(292, 2715) Crown group Skew-t$$^{\mathrm{a,b}}$$ Gamma$$^{\mathrm{a,c}}$$ Primates ST(0.878, 0.169, 2.41, 94.6) G(78.6, 77.6) Strepsirrhini ST(0.580, 0.0567, 2.12, 149) G(262, 327) Haplorrhini ST(0.705, 0.0628, 1.66, 960) G(254, 411) Anthropoidea ST(0.415, 0.0291, 0.949, 294) G(271, 363) Catarrhini ST(0.292, 0.0206, 0.995, 167) G(316, 733) Hominoidea ST(0.185, 0.0167, 2.44, 312) G(316, 1040) Human-Gorilla ST(0.0996, 0.0103, 19.6, 100) G(311, 1575) Human-Chimp ST(0.0788, 0.00687, 3.51, 6.15) G(292, 2715) $$^{\mathrm{a}}$$Densities calculated under fossil calibration strategy A and AR model using a time unit of 100 My. $$^{\mathrm{b}}$$The parameters of the skew-t distribution are location, scale, shape, and df. $$^{\mathrm{c}}$$Note that here we use the shape ($$\alpha )$$ and rate ($$\beta )$$ parameterization. For the scale parameterization use $$s=1$$/$$\beta$$. The mean is $$\alpha$$/$$\beta$$ and variance is $$\alpha$$/$$\beta^{\mathrm{2}}$$. Effect of the Clock Model The clock model has a strong impact on posterior time estimates, particularly for the most ancient nodes in the phylogeny. Under the IR model, the ages of Euarchonta, Primatomorpha, Primates, Haplorrhini, and Strepsirrhini are substantially older than those estimated under the AR model (Fig. 5B). Posterior means and 95% CIs for locus (partition) rates obtained under both clock models are given in the Supplementary Material spreadsheet available on Dryad. Results of Bayesian model selection of clock model using thermodynamic integration are shown in Table 5. The AR model has the highest marginal likelihood in 5 of the 6 partitions analysed, with the posterior model probability > 90% in two partitions, and 79%, 66%, 53%, and 29% in the other four. When the six partitions are analysed in a multipartition data set, the posterior probability is virtually 100% in favor of the AR model. We note that ideally, the marginal likelihood calculations should have been carried out on the complete data set, but unfortunately, this is so computationally expensive that it cannot be done in a feasible amount of time. Thus, further work is necessary to confirm whether the preference for the AR model will remain in analysis of a more taxonomically dense primate phylogeny. Table 5. Bayesian selection of relaxed-clock model Data set Model Log marginal $$L^{\mathrm{a}}$$ BF$$^{\mathrm{b}}$$ $$P^{\mathrm{c}}$$ Mitochondrial 1st and 2nd c.p. SC –16,519.03 (0.010) $$1.3\times 10^{\mathrm{-18}}$$ $$1.2\times 10^{\mathrm{-18}}$$ IR –16,480.58 (0.021) 0.063 0.060 AR –16,477.82 (0.035) — 0.94 Mitochondrial 3rd c.p. SC –16,684.50 (0.014) — 0.61 IR –16,686.29 (0.043) 0.17 0.10 AR –16,685.26 (0.040) 0.47 0.29 Mitochondrial RNA SC –7906.85 (0.0087) 0.74 0.39 IR –7908.40 (0.015) 0.16 0.08 AR –7906.55 (0.023) — 0.53 Nuclear 1st and 2nd c.p. SC –32,179.80 (0.0092) 0.0047 0.0037 IR –32,175.77 (0.022) 0.27 0.21 AR –32,174.44 (0.032) — 0.79 Nuclear 3rd c.p. SC –24,535.33 (0.012) $$7.2\times 10^{\mathrm{-12}}$$ $$6.7\times 10^{\mathrm{-12}}$$ IR –24,512.45 (0.038) 0.062 0.058 AR –24,509.67 (0.030) — 0.94 Nuclear UTR and introns SC –64,739.20 (0.016) $$5.7\times 10^{\mathrm{-4}}$$ $$3.8\times 10^{\mathrm{-4}}$$ IR –64,732.41 (0.038) 0.51 0.34 AR –64,731.73 (0.046) — 0.66 All$$^{\mathrm{d}}$$ SC –162,684.8 (0.024) $$2.1 \times 10^{\mathrm{-103}}$$ $$2.1\times 10^{\mathrm{-103}}$$ IR –162,467.0 (0.086) $$8.4\times 10^{\mathrm{-9}}$$ $$8.4\times 10^{\mathrm{-9}}$$ AR –162,448.4 (0.15) — 1.00 Data set Model Log marginal $$L^{\mathrm{a}}$$ BF$$^{\mathrm{b}}$$ $$P^{\mathrm{c}}$$ Mitochondrial 1st and 2nd c.p. SC –16,519.03 (0.010) $$1.3\times 10^{\mathrm{-18}}$$ $$1.2\times 10^{\mathrm{-18}}$$ IR –16,480.58 (0.021) 0.063 0.060 AR –16,477.82 (0.035) — 0.94 Mitochondrial 3rd c.p. SC –16,684.50 (0.014) — 0.61 IR –16,686.29 (0.043) 0.17 0.10 AR –16,685.26 (0.040) 0.47 0.29 Mitochondrial RNA SC –7906.85 (0.0087) 0.74 0.39 IR –7908.40 (0.015) 0.16 0.08 AR –7906.55 (0.023) — 0.53 Nuclear 1st and 2nd c.p. SC –32,179.80 (0.0092) 0.0047 0.0037 IR –32,175.77 (0.022) 0.27 0.21 AR –32,174.44 (0.032) — 0.79 Nuclear 3rd c.p. SC –24,535.33 (0.012) $$7.2\times 10^{\mathrm{-12}}$$ $$6.7\times 10^{\mathrm{-12}}$$ IR –24,512.45 (0.038) 0.062 0.058 AR –24,509.67 (0.030) — 0.94 Nuclear UTR and introns SC –64,739.20 (0.016) $$5.7\times 10^{\mathrm{-4}}$$ $$3.8\times 10^{\mathrm{-4}}$$ IR –64,732.41 (0.038) 0.51 0.34 AR –64,731.73 (0.046) — 0.66 All$$^{\mathrm{d}}$$ SC –162,684.8 (0.024) $$2.1 \times 10^{\mathrm{-103}}$$ $$2.1\times 10^{\mathrm{-103}}$$ IR –162,467.0 (0.086) $$8.4\times 10^{\mathrm{-9}}$$ $$8.4\times 10^{\mathrm{-9}}$$ AR –162,448.4 (0.15) — 1.00 Notes: Marginal likelihoods are estimated by thermodynamic integration with 64 points. The substitution model is model is HKY$$+$$G. The age of the root is fixed to one (i.e., we use a “B(0.99, 1.01)” calibration on the root in MCMCTree). The rate priors are G(2, 1) and G(2, 20) for mitochondrial and nuclear data, respectively. The prior on $$\sigma^{\mathrm{2}}$$ is G(1, 1) in all cases. The model with the highest posterior probability in each data set is shown in bold type. SC $$=$$ strict clock; IR $$=$$ independent log-normal; AR $$=$$ autocorrelated rates.$$\sigma$$ $$^{\mathrm{a}}$$Values in brackets are the standard errors (see Appendix 2). $$^{\mathrm{b}}$$The values are the BF of the given model compared with the model with highest marginal likelihood (see Appendix 2). $$^{\mathrm{c}}$$Posterior model probabilities are calculated assuming a uniform prior on models. $$^{\mathrm{d}}$$The six data sets are analysed together as six partitions. Table 5. Bayesian selection of relaxed-clock model Data set Model Log marginal $$L^{\mathrm{a}}$$ BF$$^{\mathrm{b}}$$ $$P^{\mathrm{c}}$$ Mitochondrial 1st and 2nd c.p. SC –16,519.03 (0.010) $$1.3\times 10^{\mathrm{-18}}$$ $$1.2\times 10^{\mathrm{-18}}$$ IR –16,480.58 (0.021) 0.063 0.060 AR –16,477.82 (0.035) — 0.94 Mitochondrial 3rd c.p. SC –16,684.50 (0.014) — 0.61 IR –16,686.29 (0.043) 0.17 0.10 AR –16,685.26 (0.040) 0.47 0.29 Mitochondrial RNA SC –7906.85 (0.0087) 0.74 0.39 IR –7908.40 (0.015) 0.16 0.08 AR –7906.55 (0.023) — 0.53 Nuclear 1st and 2nd c.p. SC –32,179.80 (0.0092) 0.0047 0.0037 IR –32,175.77 (0.022) 0.27 0.21 AR –32,174.44 (0.032) — 0.79 Nuclear 3rd c.p. SC –24,535.33 (0.012) $$7.2\times 10^{\mathrm{-12}}$$ $$6.7\times 10^{\mathrm{-12}}$$ IR –24,512.45 (0.038) 0.062 0.058 AR –24,509.67 (0.030) — 0.94 Nuclear UTR and introns SC –64,739.20 (0.016) $$5.7\times 10^{\mathrm{-4}}$$ $$3.8\times 10^{\mathrm{-4}}$$ IR –64,732.41 (0.038) 0.51 0.34 AR –64,731.73 (0.046) — 0.66 All$$^{\mathrm{d}}$$ SC –162,684.8 (0.024) $$2.1 \times 10^{\mathrm{-103}}$$ $$2.1\times 10^{\mathrm{-103}}$$ IR –162,467.0 (0.086) $$8.4\times 10^{\mathrm{-9}}$$ $$8.4\times 10^{\mathrm{-9}}$$ AR –162,448.4 (0.15) — 1.00 Data set Model Log marginal $$L^{\mathrm{a}}$$ BF$$^{\mathrm{b}}$$ $$P^{\mathrm{c}}$$ Mitochondrial 1st and 2nd c.p. SC –16,519.03 (0.010) $$1.3\times 10^{\mathrm{-18}}$$ $$1.2\times 10^{\mathrm{-18}}$$ IR –16,480.58 (0.021) 0.063 0.060 AR –16,477.82 (0.035) — 0.94 Mitochondrial 3rd c.p. SC –16,684.50 (0.014) — 0.61 IR –16,686.29 (0.043) 0.17 0.10 AR –16,685.26 (0.040) 0.47 0.29 Mitochondrial RNA SC –7906.85 (0.0087) 0.74 0.39 IR –7908.40 (0.015) 0.16 0.08 AR –7906.55 (0.023) — 0.53 Nuclear 1st and 2nd c.p. SC –32,179.80 (0.0092) 0.0047 0.0037 IR –32,175.77 (0.022) 0.27 0.21 AR –32,174.44 (0.032) — 0.79 Nuclear 3rd c.p. SC –24,535.33 (0.012) $$7.2\times 10^{\mathrm{-12}}$$ $$6.7\times 10^{\mathrm{-12}}$$ IR –24,512.45 (0.038) 0.062 0.058 AR –24,509.67 (0.030) — 0.94 Nuclear UTR and introns SC –64,739.20 (0.016) $$5.7\times 10^{\mathrm{-4}}$$ $$3.8\times 10^{\mathrm{-4}}$$ IR –64,732.41 (0.038) 0.51 0.34 AR –64,731.73 (0.046) — 0.66 All$$^{\mathrm{d}}$$ SC –162,684.8 (0.024) $$2.1 \times 10^{\mathrm{-103}}$$ $$2.1\times 10^{\mathrm{-103}}$$ IR –162,467.0 (0.086) $$8.4\times 10^{\mathrm{-9}}$$ $$8.4\times 10^{\mathrm{-9}}$$ AR –162,448.4 (0.15) — 1.00 Notes: Marginal likelihoods are estimated by thermodynamic integration with 64 points. The substitution model is model is HKY$$+$$G. The age of the root is fixed to one (i.e., we use a “B(0.99, 1.01)” calibration on the root in MCMCTree). The rate priors are G(2, 1) and G(2, 20) for mitochondrial and nuclear data, respectively. The prior on $$\sigma^{\mathrm{2}}$$ is G(1, 1) in all cases. The model with the highest posterior probability in each data set is shown in bold type. SC $$=$$ strict clock; IR $$=$$ independent log-normal; AR $$=$$ autocorrelated rates.$$\sigma$$ $$^{\mathrm{a}}$$Values in brackets are the standard errors (see Appendix 2). $$^{\mathrm{b}}$$The values are the BF of the given model compared with the model with highest marginal likelihood (see Appendix 2). $$^{\mathrm{c}}$$Posterior model probabilities are calculated assuming a uniform prior on models. $$^{\mathrm{d}}$$The six data sets are analysed together as six partitions. Effect of Genome-Scale Data Figure 6 shows the infinite-sites plot for the primate data analysed here. For calibration strategy A, the eight primate nodes shared between the 10-species and 372-species trees (i.e., the nodes constrained by the large genome-scale alignment, Table 2) fall in an almost perfectly straight line ($$R \quad =$$ 0.992, Fig. 6A). This indicates that for these nodes, uncertainty in the time estimates is dominated by uncertainties in the fossil calibrations rather than by uncertainties in the molecular data. For strategy A, a regression line through the origin fitted to the eight data points is $$w \quad =$$ 0.128$$t$$, meaning that for every 1 My of divergence, 0.128 My are added to the CI width (Fig. 6A). On the other hand, when considering all 371 nodes in the tree, the relationship between CI width and mean times is far from linear, and the level of uncertainty is much higher. In this case, 0.277 My are added to the CI width for every 1 My of divergence. The trend is similar under calibration strategy B (Fig. 6B), albeit in this case there is in general more uncertainty in time estimates (i.e., the slope of the regression lines is larger). This appears due to strategy B being more conservative than strategy A, that is, some of the calibration densities used in B are substantially wider, encompassing larger time ranges (Supplementary Figure S1 available on Dryad). Figure 6. View largeDownload slide Infinite-sites plot. Posterior CI width is plotted against mean posterior divergence times for A) analysis under calibration strategy A and AR clock, and B) analysis under calibration strategy B and AR clock. In both cases, black dots indicate the eight primate nodes shared between the 10-species and 372-species trees, while the grey dots represent the rest of the nodes. Solid line: regression through the origin fitted to the black dots. Dashed line: regression through the origin fitted to all the dots. Figure 6. View largeDownload slide Infinite-sites plot. Posterior CI width is plotted against mean posterior divergence times for A) analysis under calibration strategy A and AR clock, and B) analysis under calibration strategy B and AR clock. In both cases, black dots indicate the eight primate nodes shared between the 10-species and 372-species trees, while the grey dots represent the rest of the nodes. Solid line: regression through the origin fitted to the black dots. Dashed line: regression through the origin fitted to all the dots. Discussion A Phylogenomic View of Primate Divergences A primary aim of this work was to study the effect of genome-scale data on divergence time estimates on a species-level phylogeny. Given the wide availability of whole genome data for a core-set of species, it is important to know whether the use of these data for a subsample of lineages will be enough to reduce time uncertainties in a species-level phylogeny to the theoretical limit. The results of Figure 6 clearly indicate this is not the case. Although for the core ancestral nodes in the primate phylogeny, the genome-scale alignments do constrain uncertainty in time estimates close to their theoretical limit (so it is highly unlikely that adding additional molecular data for these species will improve time estimates appreciably), for species without genome-scale data, there are still substantial uncertainties left for family-level and genus-level divergences in the tree. For some nodes, the CI-width is almost as large as the node age (e.g., for Tarsiidae, the node age is 23.8 Ma with CI-width 18 My, which is 76% of the node age, see Supplementary Material spreadsheet available on Dryad). Thus, much work is still needed in order to improve time estimates for the hundreds of more recent divergences in the tree. Furthermore, application of morphological-based models for dating (Ronquist et al. 2012) and the fossilized birth–death process (Heath et al. 2014) also offer exciting prospects and challenges in obtaining time estimates for the species-level divergences (O’Reilly et al. 2015; dos Reis et al. 2016). Improving these estimates will be important in studies of primate diversification rates and to correlate primate diversification events with major geological events in the history of the Earth (such as glaciations, continental drift, the closure the Panama isthmus, etc.). Sequential Bayesian Analysis Versus Secondary Calibrations In this work, we used the posterior of times obtained under a small data set as the prior of times in a second analysis under a large data set. This approach is justified as long as the data sets are independent under the likelihood model and as long as the data sets do not overlap (i.e., they share no genes). The use of the posterior in an analysis as the prior for the next is a well-known feature of Bayesian inference (Gelman et al. 2013). Consider data that can be split into two subsets, $$D = (D_{\mathrm{1}}$$, $$D_{\mathrm{2}})$$, which are independent under the likelihood model. The posterior distribution for parameter $$\theta$$ is ${\begin{array}{*{20}l} {f(\theta \vert D)} & \propto\, {f(\theta )f(D_{1} \vert \theta )f(D_{2} \vert \theta )}\\ & \propto\, {f(\theta \vert D_{1} )f(D_{2} \vert \theta )} \hfill \\ \end{array}},$ where $$f(\theta \vert D_{1} )\propto f(\theta )f(D_{1} \vert \theta )$$ is the posterior distribution of $$\theta$$ when only $$D_{\mathrm{1}}$$ are analysed. It is apparent that using $$f(\theta \vert D_{1} )$$ as the prior when analysing $$D_{\mathrm{2}}$$ leads to the posterior for the joint data $$D$$. In other words, performing the analysis in one step (joint analysis of $$D_{\mathrm{1}}$$ and $$D_{\mathrm{2}})$$ or in two steps (posterior under $$D_{\mathrm{1}}$$ as prior under $$D_{\mathrm{2}})$$ results in the same posterior distribution. The approach we used here to analyse the primate data is justified because the likelihood model assumes that the sequence partitions are nonoverlapping and independent. However, our approach is approximate. In multiparameter models, the posterior is a multidimensional distribution that may have a complex correlation structure. Here, we ignored the correlation structure of the nine times estimated using the genomic data and approximated the corresponding high-dimensional time posterior as the product of the marginal densities of the times, with a truncation applied to ensure that descendants are younger than ancestors. Note that joint analysis of all the partitions would have been preferable, but it is computationally prohibitive. This Bayesian sequential analysis is different from the use of secondary calibrations in some dating studies (Graur and Martin 2004), where the secondary calibrations were used as point calibrations (with the uncertainties on node age estimates ignored), and where in many cases the data analysed under the secondary calibration was the same as the data analysed to obtain the calibration in the first place (Graur and Martin 2004). Clock Model An interesting result from our study is the finding that the AR model fits the primate data better than the IR relaxed-clock model. In the context of previous studies, Lepage et al. (2007) found, using BF and no fossil calibrations, that two AR relaxed clocks (CIR and log-normal) fitted real data (eukaryotes, mammals, and vertebrates) better than IR models. More recently, Lartillot et al. (2016) introduced a mixed relaxed clock that has autocorrelated- and independent-rates components. In their analysis, the mixed clock appeared to provide a better description of rate evolution in the mammal phylogeny’ however, they did not assess clock model fit with BF. Linder et al. (2011) found, also by using BF, that IR models better fit an angiosperm phylogeny better than AR models. Additionally, they found that, when analysed without fossil calibrations, the AR model fit an ape phylogeny better than the IR model. However, when analysed with fossil calibrations, the IR model fit the ape data better. In the AR model the variance of the log-rate for branches is proportional to the time of divergence, so that the variance is expected to be close to zero for closely related species. In other words, the AR model allows for “local clocks” for closely related species, while allowing the rate to drift substantially across distantly related clades. This model is, from a biological point of view, quite appealing intuitively, and it also fits anecdotal evidence where the strict clock cannot be statistically rejected among very closely related species, for example, among the apes (dos Reis et al. 2016, box 2). In contrast, the IR model assumes that the variance of the branch rates is the same for different time scales. This would appear biologically unrealistic. Arguments have been put forward in favor and against both of the two types of relaxed-clock models examined by our study (Thorne et al. 1998; Drummond et al. 2006; Ho 2009), and clearly further research is still needed to understand the merits of each clock model and how they relate to evolutionary patterns in different genomic regions (Ho 2014). This will be a challenging task given how difficult it has been to distinguish between the two models in simulation studies (Heath et al. 2012; Ho et al. 2015; Lepage et al. 2007). Five Decades of Primate Molecular Timetrees Prior to the study reported here, the estimated age of the living primate clade has spanned a 30 My differential, ranging from as young as 55 Ma (Koop et al. 1989) to as old as 87 Ma (Perelman et al. 2011). The new millennium has been a particularly active time for primate divergence time analysis. Beginning in the early 2000’s, published studies have employed a diverse assortment of data sets applied to the problem (e.g., nuclear, mitochondrial, and their combination), as well as a range of statistical methods and calibration densities. Despite this array of data and methods, all of these studies—with only one notable outlier (Chatterjee et al. 2009)—have consistently indicated that the crown Primates clade originated prior to the K–Pg event (see also Steiper and Seiffert 2012). Given the continued dearth of fossil data to support this hypothesis, however, the result continues to be viewed with scepticism by the paleoanthropological community (Bloch et al. 2007; Silcox 2008; O’Leary et al. 2013; but see Martin et al. 2007). As described at length above, the current study gives added weight to the conclusion that primates are an ancient clade of placental mammals, arising just prior to or millions of years before the K–Pg. And even though lineage diversification within the major subclades is hypothesized not to have occurred until after the commencement of the Paleogene, the separation of tarsiers from other haplorrhines, and the divergence of haplorrhines and strepsirrhines, consistently appear to proceed or nearly coincide with the K–Pg. Given that this event was unequivocally one of the most disruptive and destructive geological episodes in Earth history, the temporal coincidence speaks both to the ecological flexibility and to the evolutionary opportunism of the earliest primates. Although now extinct in North America and in Europe, the primate fossil record shows that the clade was once nearly pan-global, even potentially including Antarctica. Thus, when viewed in the context of divergence date estimates, all of which fall within a temporal window when, as now, continental and island residences would have already been sundered by significant oceanic barriers (most notably, the separation of South America from Africa by the Atlantic Ocean), we must conclude that early primates would have been able dispersers. In fact, the ability to cross barriers, both terrestrial and aquatic, and to successfully colonize new land masses, are distinct hallmarks of the primate radiation (Gingerich 1981; Yoder and Nowak 2006; de Queiroz 2005, 2014; Seiffert 2012; Beard 2016; Bloch et al. 2016). Supplementary material Data available from the Dryad Digital Repository: http://dx.doi.org/10.5061/dryad.c020q. An R package, mcmc3r, that helps the user with marginal likelihood estimation with MCMCTree is available from https://github.com/dosreislab/mcmc3r. Funding This work was supported by grant BB/J009709/1 from the Biotechnology and Biosciences Research Council (UK). Part of this work was carried out while MdR was visiting the National Evolutionary Synthesis Center (NESCent, National Science Foundation #EF-0905606) in summer 2013. joint scholarship from University College London and the Government of Mexico’s CONACYT to J.B.M.; Grant from the Burroughs Wellcome Fund and by DEB-1354610 from the National Science Foundation to A.D.Y. Dedication This article is dedicated to the memory of our co-author Gregg F. Gunnell and to the lasting value of his work. Appendix 1 Justifications for dates assigned to 17 fossil calibrations in this study (Table 3) are given below: 8 calibrations for strategy A (SA), 1 calibration for strategy B (SB), and 8 calibrations shared among both strategies (SAB). The justifications for the remaining calibrations in Table 3 are given in dos Reis et al. (2012; see also Benton et al. 2009). In some cases, the dates used are not exactly those published in cited references. In these cases, the dates utilized reflect published as well as unpublished information or adjustments deemed necessary given the uncertainty of some dates. In any case, the discrepancies are always small and unimportant considering the breadth of the fossil calibrations. Note that specifying maximum bounds is a difficult task, in particular because absence of fossil evidence is not evidence that a clade did not exist in a point in time (Ho and Philips 2009). We use stem fossils as benchmark points onto which to construct diffuse maximum bounds (i.e., with a large probability of violation, $$p_{\mathrm{U}})$$ on some node ages. Hominini | Homo-Pan | 7.5 Ma–10 Ma | SA The minimum age for the divergence of hominins is placed at 7.5 Ma and is based on the appearance of †Sahelanthropus at 7.2 Ma (Brunet et al. 2002, 2005; Lebatard et al. 2008). There is some controversy as to the proper taxonomic position of Sahelanthropus (Wolpoff et al. 2006; MacLatchy et al. 2010), but we regard it as the oldest record of a plausible crown hominin. Sahelanthropus comes from the Anthracotheriid Unit of an unnamed formation in the Mega-Chad Basin in Chad (Brunet et al. 2005). The associated mammalian fauna is very similar to that found from the Nawata Formation, Lothagam, Kenya which may be as old as 7.4 Ma (MacDougall 2003). The divergence of the hominin lineage seems unlikely to have occurred before the appearance of the potential gorillin Chororapithecus at 10 Ma (Suwa et al. 2007; Harrison 2010a). Homininae | Gorilla-Homo | 10 Ma–13.2 Ma | SA The minimum age for the divergence of crown hominines is placed at 10 Ma based on the appearance of the potential gorillin Chororapithecus. Like Sahelanthropus, the taxonomic status of Chororapithecus is not without controversy (Suwa et al. 2007; Harrison 2010a). Harrison (2010a) regards Chororapithecus as best interpreted as a stem hominin or even a stem hominid—we feel that the features that do support a relationship with gorillas are well enough established to use the date of appearance of Chororapithecus as a minimum divergence date for hominines. Chororapithecus comes from the late Miocene Beticha section of the Chorora Formation in Ethiopia and is dated at 10–10.5 Ma (Geraads et al. 2002). We use a maximum divergence date of 13.2 (Raza et al. 1983) for Sivapithecus, but it is now evident that this date might be slightly too old. In light of this, the divergence of the hominine lineage is unlikely to have taken place before the earliest appearance of the probable crown pongine Sivapithecus at 12.5 Ma (Begun 2010; Begun et al. 2012). Hominidae | Pongo-Homo | 11.2 Ma–28 Ma | SA The minimum age for the divergence of crown hominids is placed at 11.2 Ma based on the earliest appearance of the crown pongine Sivapithecus (Kappelman et al. 1991; Begun 2010; Begun et al. 2012). Sivapithecus is known from Siwalik group rocks (Chinji, Nagri, and Dhok Pathan formations) in Indo-Pakistan that range in age from 14 Ma to 5.5 Ma (Badgley and Behrensmeyer 1995) with Sivapithecus restricted to a range of 12.5 Ma to 7.4 Ma (Flynn et al. 1995). The divergence of crown hominids is unlikely to have occurred before the first appearance of Kamoyapithecus at 25 Ma (Seiffert 2010). We use 28 Ma as a slightly more conservative maximum. Catarrhini | Homo-Macaca | 25 Ma–33.7 Ma | SA The presence of the crown hominoid Kamoyapithecus (Zalmout et al. 2010) indicates that the minimum divergence time for crown Catarrhini is 25 Ma (Seiffert 2010). Kamoyapithecus is only known from the Erageliet beds, Kalakol Basalts locality of Lothidok in Kenya (Madden 1980; Leakey et al. 1995; Rasmussen and Gutierrez 2009). A soft maximum of 33.7 Ma on the age of crown catarrhines is given due to the absence of hominoids before 33.7 Ma. Anthropoidea | Catarrhini–Platyrrhini | 41 Ma–62.1 Ma | SA Calibration density constructed from fossil modeling. The effects of the K–Pg extinction are included in the model. See Wilkinson et al. (2011) for details. Haplorrhini | crown Tarsius | 45 Ma | SB The presence of the crown tarsiid Xanthorhysis in Shanxi Province, China (Beard 1998) and apparently of the genus Tarsius in fissure fills at Shanghuang in Jiangsu Province, China (Beard et al. 1994), both dating to the late middle Eocene (40–45 Ma), circumscribe the minimum divergence time for crown Haplorrhini at 45 Ma. Strepsirrhini | Lorisiformes–Lemuriformes | 37 Ma–58 MA | SA The minimum age for the divergence of crown Strepsirrhini is placed at 37 Ma based on the first appearance of the crown lorisiform Saharagalago (Seiffert et al. 2003). Saharagalago is only known from Fayum Quarry BQ-2 in the Birket Qarun Formation, Egypt. The divergence of crown strepsirrhines is unlikely to have occurred before the first appearance of the basal primate Altiatlasius from Ouarzazate in Morocco which is considered to represent the late Paleocene (Thanetian) and dating to around 58 Ma (Gheerbrant et al. 1993, 1998; Seiffert 2010). Primates | Haplorrhini–Strepsirrhini | 57.6 Ma–88.6 Ma | SA Calibration density constructed from fossil modeling. Includes the effects of the K–Pg extinction in the model. See Wilkinson et al. (2011) for details. Euarchonta | Scandentia–Primates | 65 Ma–130 Ma | SA The minimum age for the divergence of crown Euarchonta is placed at 65 Ma based on the first appearance of the crown euarchontan Purgatorius (Bloch et al. 2007). Purgatorius is known from the early Paleocene (Puercan) Tullock and Bear Formations in Montana (Clemens 1974; Buckley 1997; Clemens 2004; Chester et al. 2015) and from the earliest Paleocene Ravenscrag Formation in Saskatchewan (Fox and Scott 2011). The divergence of Euarchonta is unlikely to have been before the appearance of placental mammals by at least 130 Ma (Luo 2007; but see Luo et al. 2011 for a potential 130 My old eutherian). Lorisidae | Nycticebus–Perodicticus | 14 Ma–37 Ma | SAB The minimum age for the divergence of crown Lorisidae is placed at 14 Ma based on an undescribed genus and species from Fort Ternan in Kenya cited by Walker (1978) and Harrison (2010b). The minimum age could possibly be as old as 19 Ma if Mioeuoticus (Leakey in Bishop 1962; Walker 1978) represents a crown lorisid (Harrison 2010b). Fossil lorisids are known from the early to middle Miocene in Africa (Phillips and Walker 2000, 2002; Harrison 2010b) and from the late Miocene of Pakistan (Jacobs 1981). The divergence of lorisoids is unlikely to have occurred before the first appearance of the potential stem lorisid Karanisia at BQ-2 in Egypt (Seiffert et al. 2003). Galagidae | Galago–Euoticus | 15 Ma–37 Ma | SAB The minimum age for the divergence of crown Galagidae is placed at 15 Ma based on an undescribed genus and species from Maboko Island in Kenya cited by McCrossin (1999) and Harrison (2010b). The minimum age could possibly be as old as 19 Ma if either Progalago or Komba represent a crown galagid (MacInnes 1943; Simpson 1967; Harrison 2010b). Fossil galagids are known from the early Miocene through early Pleistocene in Africa (Phillips and Walker 2002; Harrison 2010b). The divergence of galagids is unlikely to have occurred before the first appearance of the potential stem lorisid Karanisia at BQ-2 in Egypt (Seiffert et al. 2003). Lorisiformes | Nycticebus–Perodicticus | 14 Ma–37 Ma | SAB The minimum age for the divergence of crown Lorisiformes is placed at 18 Ma based on the appearance of the potential crown lorisoid Mioeuoticus in the early Miocene of East Africa (Harrison 2010b). The divergence of lorisiforms is unlikely to have occurred before the first appearances of Karanisia at BQ-2 in Egypt (Seiffert et al. 2003). We use 38 Ma as a conservative soft maximum. Platyrrhini | Pitheciidae–Callitrichidae | 15.7 Ma–33 Ma | SAB The minimum age for the divergence of crown Platyrrhini is based on the first occurrence of the crown pitheciine Proteropithecia dated at 15.7 Ma (Kay et al. 1998; Fleagle and Tejedor 2002). The minimum age could be as much 18 Ma if either (or both) Soriacebus or Carlocebus represent crown pitheciins (Fleagle et al. 1987, 1995; Fleagle 1990; Bown and Fleagle 1993; Rosenberger 2011). All of these taxa are known from the early and middle Miocene of Argentina (Fleagle and Tejedor 2002). The divergence of platyrrhines is unlikely to have occurred before the appearance of the crown catarrhine Catopithecus (33 Ma, Fayum, Egypt) although a recently published report has claimed a 36 Ma date for a stem platyrrhine from Peru (Fleagle and Tejedor 2002). It remains unclear how this older date was derived, however, and requires further substantiation. Atelidae | Ateles-Alouatta | 12.8 Ma–18 Ma | SAB The minimum age for the divergence of crown Atelidae (as recognized by Rosenberger 2011; subfamily Atelinae of others) is based on the first appearance of the crown atelid Stirtonia (Hershkovitz 1970; Rosenberger 2011) at 12.8 Ma. Fossil atelids are known from the middle Miocene of Colombia and the Quaternary of Brazil and the Greater Antilles (MacPhee and Horovitz 2002). The divergence of atelids is unlikely to have occurred before the first appearance of the potential stem or crown atelid Soriacebus at 18 Ma (Bown and Fleagle 1993; Fleagle et al. 1995). Cebidae | Cebus–Saimiri | 12.8 Ma–18 Ma | SAB The minimum age for the divergence of crown Cebidae is based on the first appearance of the crown cebid Neosaimiri (Stirton 1951; Hartwig and Meldrum 2002) dated at 12.8 from La Venta, Colombia. There are several older potential crown cebids including Dolichocebus and Tremacebus from Argentina and Chilecebus from Chile, all dated to around 20 Ma but it remains unclear how these taxa relate to the crown group. Recently, an additional potential crown cebid has been described from Panama (Bloch et al. 2016) dated at 20.9 Ma, which if substantiated would push the potential maximum bound to at least 21 Ma. Here we have accepted the notion that the divergence of cebids is unlikely to have occurred before the first appearance of the potential stem or crown atelid Soriacebus at 18 Ma (Bown and Fleagle 1993; Fleagle et al. 1995) but acknowledge that this date could be as old as 21–22 Ma. Cercopithecidae | Papionini–Cercopithecini | 5 Ma–23 Ma | SAB The minimum age for the divergence of crown Cercopithecidae is based on the first appearance of Parapapio in the late Miocene (5 Ma) at Lothagam, Kenya. It is potentially possible that some specimens of Parapapio from Lothagam could be as old as 7.4 Ma (Jablonski and Frost 2010). The divergence of cercopithecids is unlikely to have occurred before the first appearance of the stem cercopithecoid Prohylobates. The oldest documented Prohylobates specimens are from Wadi Moghra in the Qatarra Depression, Egypt dated to approximately 19.5 Ma. However, we view Kamoyapithecus, dated at 25 Ma as a crown hominoid, which implies that at least stem cercopithecoids were in existence at that time. Given the controversial nature of our views on Kamoyapithecus, we have used a date of 23 Ma as the likely maximum divergence time for crown cercopithecids. Colobinae | Colobini–Presbytini | 9.8 Ma–23 Ma | SAB The minimum age for the divergence of crown Colobinae is based on the first appearance of Microcolobus in the middle Miocene (9.8 Ma) at Ngeringerowa, Kenya (Benefit and Pickford 1986; Jablonski and Frost 2010). The divergence of colobines is unlikely to have occurred before the first appearance of the stem cercopithecoid Prohylobates. The oldest documented Prohylobates specimens are from Wadi Moghra in the Qatarra Depression, Egypt dated to approximately 19.5 Ma. However, we view Kamoyapithecus, dated at 25 Ma as a crown hominoid, which implies that at least stem cercopithecoids were in existence at that time. Given the controversial nature of our views on Kamoyapithecus, we have used a date of 23 Ma as the likely maximum divergence time for crown colobines. Appendix 2 Here, we briefly describe our new implementation of BF calculation in MCMCTree. Our approach is the thermodynamic integration-Gaussian quadrature method implemented recently in the program BPP. The mathematical details are given in Rannala and Yang (2017). In this article, we compare different rate models, $$m_{i}$$, which differ only in the density of the branch rates while the prior on divergence times and the sequence likelihood are the same between the models. The posterior distribution of times (t) and rates (r) given the sequence data $$D$$, and given a clock model $$m_{i}$$ is thus $f({\rm {\bf t}},{\rm {\bf r}}\vert D,m_{i} )=\frac{1}{z_{i} }f({\rm {\bf t}})f({\rm {\bf r}}\vert {\rm {\bf t}},m_{i} )f(D\vert {\rm {\bf t}},{\rm {\bf r}}),$ where $z_{i} =\int {f({\rm {\bf t}})f({\rm {\bf r}}\vert {\rm {\bf t}},m_{i} )f(D\vert {\rm {\bf t}},{\rm {\bf r}})\mbox{d}{\rm {\bf t}}\mbox{d}{\rm {\bf r}}}$ is the marginal likelihood of the data for model $$m_{i}$$. Let $$m_{r}$$ be the model with highest marginal likelihood, and let BF$$_{ir} \quad =$$$$z_{i}$$/$$z_{r}$$ be the BF of model $$i$$ over model $$r$$. Then the posterior probability of model $$i$$ is $\Pr (m_{i} \vert D)=\frac{z_{i} \Pr (m_{i} )}{\sum\limits_j {z_{j} \Pr (m_{j} )} }=\frac{z_{i} /z_{r} \Pr (m_{i} )}{\sum\limits_j {z_{j} /z_{r} \Pr (m_{j} )} }=\frac{\mbox{BF}_{ir} \Pr (m_{i} )}{\sum\limits_j {\mbox{BF}_{jr} \Pr (m_{j} )} },$ where the sum is over all models being tested, and Pr($$m_{\mathrm{i}})$$ is the prior model probability. We calculate $$z_{i}$$ by sampling from the power posterior $f_{\beta } ({\rm {\bf t}}\mbox{,}{\rm {\bf r}}\vert D,m_{i} )\propto f({\rm {\bf t}})f({\rm {\bf r}}\vert {\rm {\bf t}},m_{i} )f(D\vert {\rm {\bf t}}\mbox{,}{\rm {\bf r}})^{\beta },\quad 0\leqslant \beta \leqslant 1$ for a given $${\beta}$$ value. We choose $${K}\, \beta_{j}$$ values to integrate between 0 and 1 according to the Gauss-Legendre quadrature rule. The estimate of log $$z_i$$ is then given by the quadrature formula $\log z_{i} \approx \frac{1}{2}\sum\limits_{j=1}^K {w_{j} } \bar{{\ell }}_{\beta_{j} } ,$ where $$w_{j}$$ are the Gauss–Legendre quadrature weights, and $$\bar{{\ell }}_{\beta_{j} }$$ is the average of log-likelihood values, $$\ell_{\beta_{j}}$$, sampled from the power posterior with$$\beta_{j}$$. The standard error of the estimate is given by $\mbox{SE}=\frac{1}{2}\sqrt {\sum\limits_{j=1}^K {w_{j}^{2} \mbox{Var}} \left( {\ell_{\beta_{j} } } \right)/\mbox{ESS}} ,$ where the variance is calculated over the MCMC sample of $$\ell_{\beta_{j}}$$, and ESS is the effective-sample size of $$\ell_{\beta_{j}}$$. References Alfaro M.E., Santini F., Brock C., Alamillo H., Dornburg A., Rabosky D.L., Carnevale G., Harmon L.J. 2009 . Nine exceptional radiations plus high turnover explain species diversity in jawed vertebrates. Proc. Natl. Acad. Sci. USA 106 : 13410 – 13414 . Google Scholar CrossRef Search ADS Andersen M.J., Shult H.T., Cibois A., Thibault J.C., Filardi C.E., Moyle R.G. 2015 . Rapid diversification and secondary sympatry in Australo-Pacific kingfishers (Aves: Alcedinidae: Todiramphus). R. Soc. Open Sci., 2 : 140375 . Google Scholar CrossRef Search ADS Aris-Brosou S., Yang Z. 2002 . Effects of models of rate evolution on estimation of divergence dates with special reference to the metazoan 18S ribosomal RNA phylogeny. Syst. Biol. 51 : 703 – 714 . Google Scholar CrossRef Search ADS PubMed Badgley C., Behrensmeyer A.K. 1995 . Two long geological records of continental ecosystems. Palaeogeogr. Palaeoclimatol. Palaeoecol. 115 : 1 – 11 . Google Scholar CrossRef Search ADS Baele G., Lemey P., Bedford T., Rambaut A., Suchard M.A., Alekseyenko A.V. 2012a . Improving the accuracy of demographic and molecular clock model comparison while accommodating phylogenetic uncertainty. Mol. Biol. Evol. 29 : 2157 – 2167 . Google Scholar CrossRef Search ADS Baele G., Li W.L.S., Drummond A.J., Suchard M.A., Lemey P. 2012b . Accurate model selection of relaxed molecular clocks in Bayesian phylogenetics. Mol. Biol. Evol. 30 : 239 – 243 . Google Scholar CrossRef Search ADS Beard K. 1998 . A new genus of Tarsiidae (Mammalia: Primates) from the middle Eocene of Shanxi Province, China, with notes on the historical biogeography of tarsiers. Bull. Carnegie Mus. Nat. Hist. 34 : 260 – 277 . Beard K.C. 2016 . Out of Asia: Anthropoid origins and the colonization of Africa. Annu. Rev. Anthropol., 45 : 13.1 – 13.15 . Google Scholar CrossRef Search ADS Beard K.C., Qi T., Dawson M.R., Wang B., Li C. 1994 . A diverse new primate fauna from middle Eocene fissure-fillings in southeastern China. Nature 368 : 604 – 609 . Google Scholar CrossRef Search ADS PubMed Begun D.R. 2010 . Miocene hominids and the origins of the African apes and humans. Annu. Rev. Anthropol. 39 : 67 – 84 . Google Scholar CrossRef Search ADS Begun D.R., Nargolwalla M.C., Kordos L. 2012 . European Miocene hominids and the origin of the African ape and human clade. Evol. Anthropol. 21 : 10 – 23 . Google Scholar CrossRef Search ADS PubMed Benefit B.R., Pickford M. 1986 . Miocene fossil cercopithecoids from Kenya. Am. J. Phys. Anthropol. 69 : 441 – 464 . Google Scholar CrossRef Search ADS Benton M.J., Donoghue P.C. 2007 . Paleontological evidence to date the tree of life. Mol. Biol. Evol. 24 : 26 – 53 . Google Scholar CrossRef Search ADS PubMed Benton M.J., Donoghue P.C.J., Asher R.J. 2009 . Calibrating and constraining molecular clocks. In: Hedges B.S., Kumar S., editors. The timetree of life. Oxford, England : Oxford University Press . p. 35 – 86 . Bininda-Emonds O.R., Cardillo M., Jones K.E., MacPhee R.D., Beck R.M., Grenyer R., Price S.A., Vos R.A., Gittleman J.L., Purvis A. 2007 . The delayed rise of present-day mammals. Nature 446 : 507 – 512 . Google Scholar CrossRef Search ADS PubMed Bishop W. 1962 . The mammalian fauna and geomorphological relations of the Napak volcanics, Karamoja. Geol. Surv. Uganda Rec. 1957–1958 : 1 – 18 . Bloch J.I., Silcox M.T., Boyer D.M., Sargis E.J. 2007 . New Paleocene skeletons and the relationship of plesiadapiforms to crown-clade primates. Proc. Natl. Acad. Sci. USA 104 : 1159 – 1164 . Google Scholar CrossRef Search ADS Bloch J.I., Woodruff E.D., Wood A.R., Rincon A.F., Harrington A.R., Morgan, Gary S. Foster, David A. Camilo Montes, Jaramillo, Carlos A. Jud, Nathan A. Jones Douglas S. & MacFadden Bruce J. 2016 . First North American fossil monkey and early Miocene tropical biotic interchange. Nature 533 : 243 – 246 . Google Scholar CrossRef Search ADS PubMed Bond M., Tejedor M.F., Campbell K.E. Jr., Chornogubsky L., Novo N., Goin F. 2015 . Eocene primates of South America and the African origins of New World monkeys. Nature 520 : 538 – 541 . Google Scholar CrossRef Search ADS PubMed Bown T.M., Fleagle J.G. 1993 . Systematics, biostratigraphy, and dental evolution of the Palaeothentidae, later Oligocene to early-middle Miocene (Deseadan-Santacrucian) caenolestoid marsupials of South America. J. Paleontol. 2 Suppl : 1 – 76 . Brunet M., Guy F., Pilbeam D., Lieberman D.E., Likius A., Mackaye H.T., Ponce de Leon M.S., Zollikofer C.P., Vignaud P. 2005 . New material of the earliest hominid from the Upper Miocene of Chad. Nature 434 : 752 – 755 . Google Scholar CrossRef Search ADS PubMed Brunet M., Guy F., Pilbeam D., Mackaye H.T., Likius A., Ahounta D., Beauvilain A., Blondel C., Bocherens H., Boisserie J.-R., Louis De Bonis, Yves Coppens, Jean Dejax, Christiane Denys, Philippe Duringer, Véra Eisenmann, Gongdibé Fanone, Pierre Fronty, Denis Geraads, Thomas Lehmann, Fabrice Lihoreau, Antoine Louchart, Adoum Mahamat, Gildas Merceron, Guy Mouchelin, Olga Otero, Pablo Pelaez Campomanes, Marcia Ponce De Leon, Jean-Claude Rage, Michel Sapanet, Mathieu Schuster, Jean Sudre, Pascal Tassy, Xavier Valentin, Patrick Vignaud, Laurent Viriot, Antoine Zazzo, Christoph Zollikofer. 2002 . A new hominid from the upper Miocene of Chad, Central Africa. Nature 418 : 145 – 151 . Google Scholar CrossRef Search ADS PubMed Buckley G.A. 1997 . A new species of Purgatorius (Mammalia; Primatomorpha) from the lower Paleocene Bear formation, Crazy Mountains basin, south-central Montana. J. Paleontol. 71 : 149 – 155 . Google Scholar CrossRef Search ADS Burgess R., Yang Z. 2008 . Estimation of hominoid ancestral population sizes under Bayesian coalescent models incorporating mutation rate variation and sequencing errors. Mol. Biol. Evol. 25 : 1979 – 1994 . Google Scholar CrossRef Search ADS PubMed Chatterjee H.J., Ho S.Y., Barnes I., Groves C. 2009 . Estimating the phylogeny and divergence times of primates using a supermatrix approach. BMC Evol. Biol. 9 : 259 . Google Scholar CrossRef Search ADS PubMed Chester S.G., Bloch J.I., Boyer D.M., Clemens W.A. 2015 . Oldest known euarchontan tarsals and affinities of Paleocene Purgatorius to Primates. Proc. Natl. Acad. Sci. USA 112 : 1487 – 1492 . Clemens W.A. 1974 . Purgatorius, an early paromomyid primate (Mammalia). Science 184 : 903 – 905 . Google Scholar CrossRef Search ADS PubMed Clemens W.A. 2004 . Purgatorius (Plesiadapiformes, Primates, Mammalia), a Paleocene immigrant into northeastern Montana: stratigraphic occurrences and incisor proportions. Bull. Carnegie Mus. Nat. His. 36 : 3 – 13 . Google Scholar CrossRef Search ADS Corruccini R.S., Baba M., Goodman M., Ciochon R.L., Cronin J.E. 1980 . Non-linear macromolecular evolution and the molecular clock. Evolution 34 : 1216 – 1219 . Google Scholar CrossRef Search ADS PubMed de Queiroz A. 2005 . The resurrection of oceanic dispersal in historical biogeography. Trends Ecol. Evol. 20 : 68 – 73 Google Scholar CrossRef Search ADS PubMed de Queiroz A. 2014 . The monkey’s voyage: how improbable journeys shaped the history of life. New York : Basic Books . Dornburg A, Beaulieu J.M., Oliver J.C., Near T.J. 2011 . Integrating fossil preservation biases in the selection of calibrations for molecular divergence time estimation. Syst. Biol. 60 : 519 – 527 . Google Scholar CrossRef Search ADS PubMed dos Reis M., Donoghue P.C., Yang Z. 2014a . Neither phylogenomic nor palaeontological data support a Palaeogene origin of placental mammals. Biol. Lett. 10 : 20131003 . Google Scholar CrossRef Search ADS dos Reis M., Donoghue P.C., Yang Z. 2016 . Bayesian molecular clock dating of species divergences in the genomics era. Nat. Rev. Genet. 17 : 71 – 80 . Google Scholar CrossRef Search ADS PubMed dos Reis M., Inoue J., Hasegawa M., Asher R., Donoghue P.C., Yang Z. 2012 . Phylogenomic data sets provide both precision and accuracy in estimating the timescale of placental mammal phylogeny. Proc. R. Soc. Lond. B Biol. Sci. 279 : 3491 – 3500 . Google Scholar CrossRef Search ADS dos Reis M., Thawornwattana Y., Angelis K., Telford M.J., Donoghue P.C., Yang Z. 2015 . Uncertainty in the timing of origin of animals and the limits of precision in molecular timescales. Curr. Biol. 25 : 2939 – 2950 . Google Scholar CrossRef Search ADS PubMed dos Reis M., Yang Z. 2011 . Approximate likelihood calculation for Bayesian estimation of divergence times. Mol. Biol. Evol. 28 : 2161 – 2172 . Google Scholar CrossRef Search ADS PubMed dos Reis M., Yang Z. 2013 . The unbearable uncertainty of Bayesian divergence time estimation. J. Syst. Evol. 51 : 30 – 43 . Google Scholar CrossRef Search ADS dos Reis M., Zhu T., Yang Z. 2014b . The impact of the rate prior on Bayesian estimation of divergence times with multiple loci. Syst. Biol. 63 : 555 – 565 . Google Scholar CrossRef Search ADS Drummond A., Ho S., Phillips M., Rambaut A. 2006 . Relaxed phylogenetics and dating with confidence. PLoS Biol. 4 : e88 . Google Scholar CrossRef Search ADS PubMed Duchene S., Lanfear R., Ho S.Y.W. 2014 . The impact of calibration and clock-model choice on molecular estimates of divergence times. Mol. Phylogenet. Evol. 78 : 277 – 289 . Google Scholar CrossRef Search ADS PubMed Eizirik E., Murphy W.J., Springer M.S., O’Brien S.J. 2004 . Molecular phylogeny and dating of early primate divergences. In: Ross C.F., Kay R.F., editors. Anthropoid origins. Boston, MA : Springer. p. 45 – 64 . Google Scholar CrossRef Search ADS Finstermeier K., Zinner D., Brameier M., Meyer M., Kreuz E., Hofreiter M., Roos C. 2013 . A mitogenomic phylogeny of living primates. PLoS One 8 : e69504 . Google Scholar CrossRef Search ADS PubMed Fleagle J.G. 1990 . New fossil platyrrhines from the pinturas formation. South Argent. J. Hum. Evol. 19 : 61 – 85 . Google Scholar CrossRef Search ADS Fleagle J.G., Bown T., Swisher C., Buckley G. 1995 . Age of the Pinturas and Santa Cruz formations. VI Congreso Argentino de Paleontología y Bioestratigrafía. Trelew, Chubut, Argentine : Museo Paleontológico Eugidio Feruglio p. 129 – 135 . Fleagle J.G., Powers D.W., Conroy G.C., Watters J.P. 1987 . New fossil platyrrhines from Santa-Cruz Province, Argentina. Folia Primatol. 48 : 65 – 77 . Google Scholar CrossRef Search ADS Fleagle J.G., Tejedor M.F. 2002 . “Early platyrrhines of southern South America” in Hartwig, W.C. (Ed). The Primate fossil record. Cambridge : Cambridge University Press, pp. 161 – 174 . Flynn L.J., Barry J.C., Morgan M.E., Pilbeam D., Jacobs L.L., Lindsay E.H. 1995 . Neogene Siwalik mammalian lineages—species longevities, rates of change, and modes of speciation. Palaeogeogr. Palaeocl. 115 : 249 – 264 . Google Scholar CrossRef Search ADS Fox R.C., Scott C.S. 2011 . A new, early Puercan (Earliest Paleocene) species of Purgatorius (Plesiadapiformes, Primates) from Saskatchewan, Canada. J. Paleontol. 85 : 537 – 548 . Google Scholar CrossRef Search ADS Gelman A., Carlin J.B., Stern H.S., Dunson D.B., Vehtari A., Rubin D.B. 2013 . Bayesian data analysis, Vol. 21, 3rd ed. Boca Raton, FL : CRC Press . Geraads D., Alemseged Z., Bellon H. 2002 . The late Miocene mammalian fauna of Chorora, Awash basin, Ethiopia: systematics, biochronology and 40K-40Ar ages of the associated volcanics. Tertiary Res. 113 – 122 . Gheerbrant E., Cappetta H., Feist M., Jaeger J-J., Sudre J., Vianey-Liaud M., Sigé B. 1993 . La succession des faunes de vertébrés d’âge paléocène supérieur et éocène inférieur dans le bassin d’Ouarzazate, Maroc. Contexte géologique, portée biostratigraphique et paléogéographique. Newsl Stratigr 33 – 58 . Gheerbrant E., Sudre J., Sen S., Abrial C., Marandat B., Sigé B., Vianey-Liaud M. 1998 . Nouvelles données sur les mammifères du Thanétien et de l’Yprésien du Bassin d’Ouarzazate (Maroc) et leur contexte stratigraphique. Palaeovertebrata 27 : 155 – 202 . Gingerich P.D. 1981 . Eocene Adapidae, paleobiogeography, and the origin of South American Platyrrhini. In: Ciochon R.L., Chiarelli A.B., editors. Evolutionary biology of the new world monkeys and continental drift. New York : Plenum Publishing Corporation. pp. 123 – 138 . Google Scholar CrossRef Search ADS Goodman M., Braunitzer G., Stangl A., Schrank B. 1983 . Evidence on human origins from haemoglobins of African apes. Nature 303 : 546 – 548 . Google Scholar CrossRef Search ADS PubMed Goodman M., Porter C.A., Czelusniak J., Page S.L., Schneider H., Shoshani J., Gunnell G., Groves C.P. 1998 . Toward a phylogenetic classification of Primates based on DNA evidence complemented by fossil evidence. Mol. Phylogenet. Evol. 9 : 585 – 598 . Google Scholar CrossRef Search ADS PubMed Graur D., Martin W. 2004 . Reading the entrails of chickens: molecular timescales of evolution and the illusion of precision. Trends Genet. 20 : 80 – 86 . Google Scholar CrossRef Search ADS PubMed Harrison T. 2010a . Dendropithecoidea, proconsuloidea, and hominoidea. In: Werdelin L., Sanders W., editors. Cenozoic mammals of Africa. Berkeley, CA : University of California Press. pp. 429 – 469 . Google Scholar CrossRef Search ADS Harrison T. 2010b . Later Tertiary lorisiformes. In: Werdelin L., Sanders W., editors. Cenozoic mammals of Africa. Berkeley, CA : University of California Press . pp. 333 – 349 . Google Scholar CrossRef Search ADS Hartwig W.C., Meldrum D.J. 2002 . Miocene platyrrhines of the northern Neotropics. Cam. S. Biol. Evol. Anthropol. 175 – 188 . Hasegawa M., Kishino H., Yano T. 1985 . Dating the human-ape splitting by a molecular clock of mitochondrial DNA. J. Mol. Evol. 22 : 160 – 174 . Google Scholar CrossRef Search ADS PubMed Heath T.A., Holder M.T., Huelsenbeck J.P. 2012 . A Dirichlet process prior for estimating lineage-specific substitution rates. Mol. Biol. Evol. 29 : 939 – 955 . Google Scholar CrossRef Search ADS PubMed Heath T.A., Huelsenbeck J.P., Stadler T. 2014 . The fossilized birth–death process for coherent calibration of divergence-time estimates. Proc. Natl. Acad. Sci. USA 111 : E2957 – E2966 . Google Scholar CrossRef Search ADS Hershkovitz P. 1970 . Notes on tertiary platyrrhine monkeys and description of a new genus from the late Miocene of Colombia. Folia Primatol. 12 : 1 – 37 . Google Scholar CrossRef Search ADS PubMed Ho S.Y. 2009 . An examination of phylogenetic models of substitution rate variation among lineages. Biol. Lett. 5 : 421 – 424 . Google Scholar CrossRef Search ADS PubMed Ho S.Y., Duchene S. 2014 . Molecular-clock methods for estimating evolutionary rates and timescales. Mol. Ecol., 23 : 5947 – 5965 . Google Scholar CrossRef Search ADS Ho S.Y., Duchene S., Duchene D. 2015 . Simulating and detecting autocorrelation of molecular evolutionary rates among lineages. Mol. Ecol. Resour. 15 : 688 – 696 . Google Scholar CrossRef Search ADS PubMed Ho S.Y., Philips M.J. 2009 . Accounting for calibration uncertainty in phylogenetic estimation of evolutionary divergence times. Syst. Biol. 58 : 367 – 380 . Google Scholar CrossRef Search ADS PubMed Ho S.Y., Phillips M.J., Drummond A.J., Cooper A. 2005 . Accuracy of rate estimation using relaxed-clock models with a critical focus on the early metazoan radiation. Mol. Biol. Evol. 22 : 1355 – 1363 . Google Scholar CrossRef Search ADS PubMed Ho, S.Y., 2014 . The changing face of the molecular clock. Trends Ecol Evol. 29 : 496 – 503 . Google Scholar CrossRef Search ADS PubMed Inoue J., Donoghue P.C.H., Yang Z. 2010 . The impact of the representation of fossil calibrations on Bayesian estimation of species divergence times. Syst. Biol. 59 : 74 – 89 . Google Scholar CrossRef Search ADS PubMed Jablonski N., Frost S. 2010 . Cercopithecoidea. In: Werdelin L., Sanders. W., editors. Cenozoic mammals of Africa. Berkeley : University of California Press . p. 393 – 428 . Google Scholar CrossRef Search ADS Jacobs L.L. 1981 . Miocene lorisid primates from the Pakistan Siwaliks. Nature 289 : 585 – 587 . Google Scholar CrossRef Search ADS PubMed Kappelman J., Kelley J., Pilbeam D., Sheikh K.A., Ward S., Anwar M., Barry J.C., Brown B., Hake P., Johnson N.M. et al. 1991 . The earliest occurrence of Sivapithecus from the middle Miocene Chinji Formation of Pakistan. J. Hum. Evol. 21 : 61 – 73 . Google Scholar CrossRef Search ADS Kay R.F., Johnson D., Meldrum D.J. 1998 . A new pitheciin primate from the middle Miocene of Argentina. Am. J. Primatol., 45 : 317 – 336 . Google Scholar CrossRef Search ADS Kishino H., Thorne J.L., Bruno W.J. 2001 . Performance of a divergence time estimation method under a probabilistic model of rate evolution. Mol. Biol. Evol. 18 : 352 – 361 . Google Scholar CrossRef Search ADS PubMed Koop B.F., Tagle D.A., Goodman M., Slightom J.L. 1989 . A molecular view of primate phylogeny and important systematic and evolutionary questions. Mol. Biol. Evol. 6 : 580 – 612 . Google Scholar PubMed Kumar S., Hedges S.B. 1998 . A molecular timescale for vertebrate evolution. Nature 392 : 917 – 920 . Google Scholar CrossRef Search ADS PubMed Lartillot N., Philippe H. 2006 . Computing Bayes factors using thermodynamic integration. Syst. Biol. 55 : 195 – 207 . Google Scholar CrossRef Search ADS PubMed Lartillot N., Phillips M.J., Ronquist F. 2016 . A mixed relaxed clock model. Phil. Trans. R. Soc. Lond. B 371 : 20150132 . Google Scholar CrossRef Search ADS Leakey M.G., Ungar P.S., Walker A. 1995 . A new genus of large primate from the late Oligocene of Lothidok, Turkana District, Kenya. J. Hum. Evol. 28 : 519 – 531 . Google Scholar CrossRef Search ADS Lebatard A.E., Bourles D.L., Duringer P., Jolivet M., Braucher R., Carcaillet J., Schuster M., Arnaud N., Monie P., Lihoreau F. et al. 2008 . Cosmogenic nuclide dating of Sahelanthropus tchadensis and Australopithecus bahrelghazali: Mio-Pliocene hominids from Chad. Proc. Natl. Acad. Sci. USA 105 : 3226 – 3231 . Google Scholar CrossRef Search ADS Lepage T., Bryant D., Philippe H., Lartillot N. 2007 . A general comparison of relaxed molecular clock models. Mol. Biol. Evol. 24 : 2669 – 2680 . Google Scholar CrossRef Search ADS PubMed Linder H.P., Hardy C.R., Rutschmann F. 2005 . Taxon sampling effects in molecular clock dating: an example from the African Restionaceae. Mol. Phylogenet. Evol. 35 : 569 – 582 . Google Scholar CrossRef Search ADS PubMed Linder M., Britton T., Sennblad B. 2011 . Evaluation of Bayesian models of substitution rate evolution: parental guidance versus mutual independence. Syst. Biol. 60 : 329 – 342 . Google Scholar CrossRef Search ADS PubMed Liu Y., Medina R., Goffinet B. 2014 . 350 my of mitochondrial genome stasis in mosses, an early land plant lineage. Mol. Biol. Evol. 31 : 2586 – 2591 . Google Scholar CrossRef Search ADS PubMed Lovejoy C.O., Burstein H., Heiple K.H. 1972 . Primate phylogeny and immunological distance. Science 176 : 803 – 805 . Google Scholar CrossRef Search ADS PubMed Luo Z.X. 2007 . Transformation and diversification in early mammal evolution. Nature 450 : 1011 – 1019 . Google Scholar CrossRef Search ADS PubMed Luo Z.X., Yuan C.X., Meng Q.J., Ji Q. 2011 . A Jurassic eutherian mammal and divergence of marsupials and placentals. Nature 476 : 442 – 445 . Google Scholar CrossRef Search ADS PubMed McDougall I.A.N., Feibel C.S. 2003 . Numerical age control for the Miocene-Pliocene succession at Lothagam, a hominoid-bearing sequence in the northern Kenya Rift. In: Leakey M.G.H., Harris J.M., editors. Lothagam: the dawn of humanity in Eastern Africa. New York : Columbia University Press. p. 43 – 64 . MacInnes D.G. 1943 . Notes on the East African primates. J. East Afr. Uganda Nat. Hist. Soc. 39 : 521 – 530 . MacLatchy L., DeSilva J., Sanders W., Wood B. 2010 . Hominini. In: Werdelin L., Sanders W. editors. Cenozoic mammals of Africa. Berkeley : University of California Press . p. 471 – 540 . Google Scholar CrossRef Search ADS MacPhee R.D., Horovitz I. 2002 . Extinct quaternary platyrrhines of the Greater Antilles and Brazil. Cam. S. Biol. Evol. Anthropol. 189 – 200 . Madden C.T. 1980 . NewProconsul (Xenopithecus) from the Miocene of Kenya. Primates 21 : 241 – 252 . Google Scholar CrossRef Search ADS Marshall C.R. 1990 . The fossil record and estimating divergence times between lineages: maximum divergence times and the importance of reliable phylogenies. J. Mol. Evol. 30 : 400 – 408 . Google Scholar CrossRef Search ADS PubMed Marshall C.R. 2008 . A simple method for bracketing absolute divergence times on molecular phylogenies using multiple fossil calibration points. Am. Nat. 171 : 726 – 742 . Google Scholar CrossRef Search ADS PubMed Martin R.D. 1993 . Primate origins: plugging the gaps. Nature 363 : 223 – 234 . Google Scholar CrossRef Search ADS PubMed Martin R.D., Soligo C., Tavaré S. 2007 . Primate origins: implications of a Cretaceous ancestry. Folia Primatol. 78 : 277 – 296 . Google Scholar CrossRef Search ADS PubMed McCrossin M.L. 1999 . Phylogenetic relationships and paleoecological adaptations of a new bushbaby from the middle Miocene of Kenya. Am. J. Phys. Anthropol. 89 : 195 – 196 . Nascimento F.F., Lazar A., SeuÃ¡nez H.N., Bovincino C.R. 2014 . Reanalysis of the biogeographical hypothesis of range expansion between robust and gracile capuchin monkeys. J. Phylogeogr. 42 : 1349 – 1357 . Nascimento F.F., dos Reis M., Yang Z. 2017 . A biologist’s guide to Bayesian phylogenetic analysis. Nat. Ecol. Evol. 1 : 1446 – 1454 . Google Scholar CrossRef Search ADS PubMed O’Leary M.A., Bloch J.I., Flynn J.J., Gaudin T.J., Giallombardo A., Norberto P. Giannini, Suzann L. Goldberg, Kraatz, Brian P. Zhe-Xi Luo, Jin Meng, Xijun Ni, Novacek, Michael J. Perini, Fernando A. Randall, Zachary S. Rougier, Guillermo W. Sargis, Eric J. Silcox, Mary T. Simmons, Nancy B. Michelle Spaulding, Velazco, Paúl M. Marcelo Weksler, Wible, John R. Cirranello. Andrea L. 2013 . The placental mammal ancestor and the post-K-Pg radiation of placentals. Science 339 : 662 – 667 . Google Scholar CrossRef Search ADS PubMed O’Reilly J., dos Reis M., Donoghue P.C.J. 2015 . Dating tips for divergence time estimation. Trends Genet. 31 : 637 – 650 . Google Scholar CrossRef Search ADS PubMed Paradis E. 2013 . Molecular dating of phylogenies by likelihood methods: a comparison of models and a new information criterion. Mol. Phylogenet. Evol. 67 : 436 – 444 . Google Scholar CrossRef Search ADS PubMed Perelman P., Johnson W.E., Roos C., SeuÃ¡nez H.N., Horvath J.E., Moreira M.A.M., Kessing B., Pontius J., Roelke M., Rumpler Y., Schneider, Maria Paula C. Artur Silva, O’Brien, Stephen J. Pecon-Slattery. Jill 2011 . A molecular phylogeny of living primates. PLoS Genet. 7 : e1001342 EP . Google Scholar CrossRef Search ADS PubMed Phillips E.M., Walker A. 2000 . A new species of fossil lorisid from the miocene of east Africa. Primates 41 : 367 – 372 . Google Scholar CrossRef Search ADS Phillips E.M., Walker A. 2002 . Fossil lorisoids. In: Hartwig W.C., editor. The primate fossil record. Cambridge : Cambridge University Press. p. 83 – 95 . Poux C., Douzery E.J. 2004 . Primate phylogeny, evolutionary rate variations, and divergence times: a contribution from the nuclear gene IRBP. Am. J. Phys. Anthropol. 124 : 1 – 16 . Google Scholar CrossRef Search ADS PubMed Pozzi L., Hodgson J.A., Burrell A.S., Sterner K.N., Raaum R.L., Disotell T.R. 2014 . Primate phylogenetic relationships and divergence dates inferred from complete mitochondrial genomes. Mol. Phylogenet. Evol. 75 : 165 – 183 . Google Scholar CrossRef Search ADS PubMed Prum R.O., Berv J.S., Dornburg A., Field D.J., Townsend J.P., Lemmon E.M., Lemmon A.R. 2015 . A comprehensive phylogeny of birds (Aves) using targeted next-generation DNA sequencing. Nature 526 : 569 – 573 . Google Scholar CrossRef Search ADS PubMed Purvis A. 1995 . A composite estimate of primate phylogeny. Phil. Trans. R. Soc. Lond. B , 348 : 405 – 421 . Google Scholar CrossRef Search ADS Radinsky L. 1978 . Do albumin clocks run on time? Science 200 : 1182 – 1183 . Google Scholar CrossRef Search ADS PubMed Rannala B. 2016 . Conceptual issues in Bayesian divergence time estimation. Phil. Trans. R. Soc. B, 371 : 20150134 . Google Scholar CrossRef Search ADS Rannala B., Yang Z. 2007 . Inferring speciation times under an episodic molecular clock. Syst. Biol. 56 : 453 – 466 . Google Scholar CrossRef Search ADS PubMed Rannala B., Yang Z. 2017 . Efficient Bayesian species tree inference under the multispecies coalescent. Syst. Biol. 66 : 823 – 842 . Google Scholar PubMed Rasmussen D.T., Gutierrez M. 2009 . A mammalian fauna from the late oligocene of Northwestern Kenya. Palaeontogr. Abt. A 288 : 1 – 52 . Raza S.M., Barry J.C., Pilbeam D., Rose M.D., Shah S.M.I., Ward S. 1983 . New hominoid primates from the middle Miocene Chinji Formation, Potwar Plateau, Pakistan. Nature 306 : 52 – 54 . Google Scholar CrossRef Search ADS Read D.W., Lestrel P.E. 1970 . Hominid phylogeny and immunology: a critical appraisal. Science 168 : 578 – 580 . Google Scholar CrossRef Search ADS PubMed Rieux A., Eriksson A., Li M., Sobkowiak B., Weinert L.A., Warmuth V., Ruiz-Linares A., Manica A., Balloux F. 2014 . Improved calibration of the human mitochondrial clock using ancient genomes. Mol. Biol. Evol. 31 : 2780 – 2792 . Google Scholar CrossRef Search ADS PubMed Ronquist F., Klopfstein S., Vilhelmsen L., Schulmeister S., Murray D.L., Rasnitsyn P. 2012 . A total-evidence approach to dating with fossils, applied to the early radiation of the Hymenoptera. Syst. Biol. 61 : 973 – 999 . Google Scholar CrossRef Search ADS PubMed Rosenberger A. 2011 . Evolutionary morphology, platyrrhine evolution, and systematics. Anat. Rec. 294 : 1955 – 1974 . Google Scholar CrossRef Search ADS Rutschmann F., Eriksson T., Abu Salim K., Conti E. 2007 . Assessing calibration uncertainty in molecular dating: the assignment of fossils to alternative calibration points. Syst. Biol. 56 : 591 – 608 . Google Scholar CrossRef Search ADS PubMed Sanderson M.J. 1997 . A nonparametric approach to estimating divergence times in the absence of rate constancy. Mol. Biol. Evol. 14 : 1218 – 1231 . Google Scholar CrossRef Search ADS Sarich V.M., Wilson A.C. 1967 . Immunological time scale for Hominoid evolution. Science 158 : 1200 – 1203 Google Scholar CrossRef Search ADS PubMed Seiffert E.R. 2010 . Chronology of Paleogene mammal localities. In: Werdelin L., Sanders W., editors. Cenozoic mammals of Africa. Berkeley : University of California Press. p. 19 – 26 . Google Scholar CrossRef Search ADS Seiffert E.R. 2012 . Early primate evolution in Afro-Arabia. Evol. Anthropol. 21 : 239 – 253 . Google Scholar CrossRef Search ADS PubMed Seiffert E.R., Simons E.L., Attia Y. 2003 . Fossil evidence for an ancient divergence of lorises and galagos. Nature 422 : 421 – 424 . Google Scholar CrossRef Search ADS PubMed Shapiro B., Rambaut A., Drummond A.J. 2006 . Choosing appropriate substitution models for the phylogenetic analysis of protein-coding sequences. Mol. Biol. Evol. 23 : 7 – 9 . Google Scholar CrossRef Search ADS PubMed Simpson G.G. 1967 . The tertiary lorisiform primates of Africa. Bull. Mus. Comp. Zool. 136 : 39 – 62 . Silcox M.T. 2008 . The biogeographic origins of primates and euprimates: East, West, North, or South of Eden? In: Dagosto M., Sargis E.J., editors. Mammalian evolutionary morphology: a tribute to Frederick S. Szalay. New York : Springer . p. 199 – 231 . Google Scholar CrossRef Search ADS Soltis P.S., Soltis D.E., Savolainen V., Crane P.R., Barraclough T.G. 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 . Google Scholar CrossRef Search ADS Springer M.S., Meredith R.W., Gatesy J., Emerling C.A., Park J., Rabosky D.L., Stadler T., Steiner C., Ryder O.A, Janečka J.E, Fisher, Colleen A. Murphy. William J. 2012 . Macroevolutionary dynamics and historical biogeography of primate diversification inferred from a species supermatrix. PLoS One 7 : e49521 . Google Scholar CrossRef Search ADS PubMed Stamatakis A. 2014 . RAxML version 8: a tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics 30 : 1312 – 1313 . Google Scholar CrossRef Search ADS PubMed Steiper M.E., Young N.M. 2006 . Primate molecular divergence dates. Mol. Phylogenet. Evol. 41 : 384 – 394 . Google Scholar CrossRef Search ADS PubMed Steiper M.E., Seiffert E.R. 2012 . Evidence for a convergent slowdown in primate molecular rates and its implications for the timing of early primate evolution. Proc. Nat. Acad. Sci. USA 109 : 6006 – 6011 . Google Scholar CrossRef Search ADS Stirton R.A. 1951 . Ceboid monkeys from the Miocene of Colombia. Bull. Univ. Calif. Pub. Geol. Sci. 28 : 315 – 356 . Suwa G., Kono R.T., Katoh S., Asfaw B., Beyene Y. 2007 . A new species of great ape from the late Miocene epoch in Ethiopia. Nature 448 : 921 – 924 . Google Scholar CrossRef Search ADS PubMed Thorne J.L., Kishino H. 2002 . Divergence time and evolutionary rate estimation with multilocus data. Syst. Biol. 51 : 689 – 702 . Google Scholar CrossRef Search ADS PubMed Thorne J.L., Kishino H., Painter I.S. 1998 . Estimating the rate of evolution of the rate of molecular evolution. Mol. Biol. Evol. 15 : 1647 – 1657 . Google Scholar CrossRef Search ADS PubMed Uzzell T., Pilbeam D. 1971 . Phyletic divergence dates of hominoid primates: a comparison of fossil and molecular data. Evolution 25 : 615 – 635 . Google Scholar CrossRef Search ADS PubMed Walker A.C. 1978 . Prosimian primates. In: Maglio V.J., Cooke H.B.S., editors. Evolution of African mammals. Cambridge : Harvard University Press. p. 90 – 99 . Warnock R.C., Parham J.F., Joyce W.G., Lyson T.R., Donoghue P.C. 2015 . Calibration uncertainty in molecular dating analyses: there is no substitute for the prior evaluation of time priors. Proc. R. Soc. B 282 : 20141013 . Google Scholar CrossRef Search ADS Wilkinson R.D., Steiper M.E., Soligo C., Martin R.D., Yang Z., Tavare S. 2011 . Dating primate divergences through an integrated analysis of palaeontological and molecular data. Syst. Biol. 60 : 16 – 31 . Google Scholar CrossRef Search ADS PubMed Wolpoff M.H., Hawks J., Senut B., Pickford M., Ahern J. 2006 . An ape or the ape: is the Toumaï cranium TM 266 a hominid. Paleoanthropology 2006 : 36 – 50 . Xie W., Lewis P.O., Fan Y., Kuo L., Chen M.-H. 2011 . Improving marginal likelihood estimation for Bayesian phylogenetic model selection. Syst. Biol. 60 : 150 – 160 . Google Scholar CrossRef Search ADS PubMed Yang Z. 1994a . Estimating the pattern of nucleotide substitution. J. Mol. Evol. 39 : 105 – 111 . Yang Z. 1994b . Maximum likelihood phylogenetic estimation from DNA sequences with variable rates over sites: approximate methods. J. Mol. Evol. 39 : 306 – 314 . Google Scholar CrossRef Search ADS Yang Z. 2007 . PAML 4: phylogenetic analysis by maximum likelihood. Mol. Biol. Evol. 24 : 1586 – 1591 . Google Scholar CrossRef Search ADS PubMed Yang Z., Rannala B. 1997 . Bayesian phylogenetic inference using DNA sequences: a Markov chain Monte Carlo Method. Mol. Biol. Evol. 14 : 717 – 724 . Google Scholar CrossRef Search ADS PubMed Yang Z., Rannala B. 2006 . Bayesian estimation of species divergence times under a molecular clock using multiple fossil calibrations with soft bounds. Mol. Biol. Evol. 23 : 212 – 226 . Google Scholar CrossRef Search ADS PubMed Yang Z., Yoder A.D. 2003 . Comparison of likelihood and Bayesian methods for estimating divergence times using multiple gene loci and calibration points, with application to a radiation of cute-looking mouse lemur species. Syst. Biol. 52 : 705 – 716 . Google Scholar CrossRef Search ADS PubMed Yoder A.D., Nowak M.D. 2006 . Has vicariance or dispersal been the predominant biogeographic force in Madagascar? Only time will tell. Annu. Rev. Ecol. Evol. Syst. 37 : 404 – 431 . Google Scholar CrossRef Search ADS Yoder A.D., Yang Z. 2004 . Divergence dates for Malagasy lemurs estimated from multiple gene loci: geological and evolutionary context. Mol. Ecol. 13 : 757 – 773 . Google Scholar CrossRef Search ADS PubMed Zalmout I.S., Sanders W.J., Maclatchy L.M., Gunnell G.F., Al-Mufarreh Y.A., Ali M.A., Nasser A.A., Al-Masari A.M., Al-Sobhi S.A., Nadhra A.O, Matari, Adel H. Wilson Jeffrey A. & Gingerich. Philip D. 2010 . New Oligocene primate from Saudi Arabia and the divergence of apes and Old World monkeys. Nature 466 : 360 – 364 . Google Scholar CrossRef Search ADS PubMed Zhang C., Stadler T., Klopfstein S., Heath T.A., Ronquist F. 2016 . Total-evidence dating under the fossilized birth-death process. Syst. Biol. 65 : 228 – 249 . Google Scholar CrossRef Search ADS PubMed Zhou X., Xu S., Xu J., Chen B., Zhou K., Yang G. 2012 . Phylogenomic analysis resolves the interordinal relationships and rapid diversification of the laurasiatherian mammals. Syst. Biol. 61 : 150 – 164 . Google Scholar CrossRef Search ADS PubMed Zhu T, dos Reis M, Yang Z. 2015 . Characterization of the uncertainty of divergence time estimation under relaxed molecular clock models using multiple loci. Syst Biol 64 : 267 – 280 . Google Scholar CrossRef Search ADS PubMed Zuckerkandl E., Pauling L. 1965 . Molecules as documents of evolutionary history. J. Theor. Biol. 8 : 357 – 366 . Google Scholar CrossRef Search ADS PubMed © The Author(s) 2018. Published by Oxford University Press, on behalf of the Society of Systematic Biologists. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited. For Permissions, please email: journals.permissions@oup.com

### Journal

Systematic BiologyOxford University Press

Published: Jan 12, 2018

## You’re reading a free preview. Subscribe to read the entire article.

### DeepDyve is your personal research library

It’s your single place to instantly
that matters to you.

over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month ### Explore the DeepDyve Library ### Search Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly ### Organize Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place. ### Access Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals. ### Your journals are on DeepDyve Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more. All the latest content is available, no embargo periods. DeepDyve ### Freelancer DeepDyve ### Pro Price FREE$49/month
\$360/year

Save searches from
PubMed

Create lists to

Export lists, citations