Royal Society Publishing


The tribe Bovini contains a number of commercially and culturally important species, such as cattle. Understanding their evolutionary time scale is important for distinguishing between post-glacial and domestication-associated population expansions, but estimates of bovine divergence times have been hindered by a lack of reliable calibration points. We present a Bayesian phylogenetic analysis of 481 mitochondrial D-loop sequences, including 228 radiocarbon-dated ancient DNA sequences, using a multi-demographic coalescent model. By employing the radiocarbon dates as internal calibrations, we co-estimate the bovine phylogeny and divergence times in a relaxed-clock framework. The analysis yields evidence for significant population expansions in both taurine and zebu cattle, European aurochs and yak clades. The divergence age estimates support domestication-associated expansion times (less than 12 kyr) for the major haplogroups of cattle. We compare the molecular and palaeontological estimates for the BisonBos divergence.

1. Introduction

Domestic animals and crops not only possess economic, cultural and religious values, but they also represent intriguing case studies for evolutionary biologists due to their distinctive selective and demographic histories. The tribe Bovini (Family Bovidae, Subfamily Bovinae) contains several species of interest, including domesticated cattle (both taurine and zebu) and yak. Abundant sequence data have been produced by recent studies of bovine domestication (e.g. Loftus et al. 1994; Bailey et al. 1996; Bradley et al. 1996; Troy et al. 2001; Bollongino et al. 2006; Edwards et al. 2007), making it possible to use sophisticated phylogenetic methods to test hypotheses about the domestication process (e.g. Finlay et al. 2007).

An accurate time scale for bovine evolution is important for distinguishing between the genetic signatures left by domestication (approx. 10 kyr ago; Dobney & Larson 2006) and post-glacial expansion following the Last Glacial Maximum (approx. 20 kyr ago; Hewitt 2000). There have been substantial disparities among the published estimates of bovine divergence times, however, leaving the time scale relatively unclear. Undoubtedly, some of these variations have arisen from differences in the phylogenetic methods, calibration points and sequence alignments that have been used (figure 1).

Figure 1

The bovine D-loop, with positions labelled according to the Bos taurus reference sequence (GenBank accession number V00654). The extended termination-associated sequence (ETAS), central domain and conserved sequence block (CSB) are labelled. Dark grey bars represent the segments analysed in the previous studies of bovines, including the 370 bp hypervariable region targeted by Bradley et al. (1996).

Overall, the bovine species complex presents a number of methodological difficulties, including the absence of reliable fossil calibrations for estimating divergence times. Deep palaeontological calibration points are possibly inappropriate for studying the evolution of such closely related species (Ho et al. 2005), while previous studies have frequently adopted published estimates for divergence times or substitution rates without due appreciation of the associated estimation error. For example, Finlay et al. (2007) used a Bayesian demographic analysis to find evidence of population expansion in domesticated bovine species but were unable to attach a reliable time scale. Although their analyses suggested that population growth occurred within the past 10 kyr, they assumed an errorless, tree-wide substitution rate of 32% per Myr, even though the original estimate came with a 95% highest posterior density (HPD) of 23–41% per Myr (Shapiro et al. 2004).

The recent proliferation of ancient DNA sequences from multiple bovine species offers an unprecedented opportunity to analyse the molecular evolutionary process in detail. Radiocarbon-dated samples can act as internal time calibrations for estimating rates and divergence times. In this study, we employ a Bayesian phylogenetic method that assigns a separate demographic model to different clades in the tree, allowing multiple species to be studied simultaneously. This method is used to analyse a dataset containing 481 mitochondrial D-loop sequences from five bovine species, including 228 radiocarbon-dated ancient DNA sequences. We produce an internally calibrated estimate of the evolutionary time scale and bovine phylogeny.

2. Material and methods

Published nucleotide sequences were obtained from the D-loop of five bovine species and were aligned manually. The dataset comprised 476 aligned sites from 481 individuals (table 1).

View this table:
Table 1

Summary of the mitochondrial D-loop alignment analysed in this study.

Bayesian phylogenetic analysis was performed using Beast v. 1.4.6 (available from using TIM+I+G substitution model as selected by comparison of Akaike Information Criterion scores. The analysis was repeated using the TIM+G model to investigate the influence of the parameter describing the proportion of invariant sites. Temporal information was incorporated into the analysis by including radiocarbon dates for ancient sequences, allowing the estimation of the age of the most recent common ancestor (MRCA) for bovine clades. A separate demographic model of exponential growth was applied, in the form of a coalescent prior, to each of the five species; the population size and growth rate were estimated in each case. An uninformative tree prior was used for the basal branches connecting the five species. An uncorrelated lognormal relaxed clock was used to accommodate possible rate variations among lineages (Drummond et al. 2006). In the Markov chain Monte Carlo analysis, samples from the posterior were drawn at every 5000 steps over a total of 50 000 000 steps, following a discarded burn-in of 5 000 000 steps. The analysis was repeated and the samples from the two runs were combined. Sampling, mixing and convergence to the stationary distribution were checked. The Beast input file, which includes sample names and ages, is available as electronic supplementary material.

3. Results

Bayesian phylogenetic analysis of the complete dataset produced a tree with 100% posterior support for the monophyly of each species (table 2), and with strong support for higher level nodes (figure 2). The inferred phylogeny confirms the results of previous mitochondrial studies, including the pairing of American bison (Bison bison) with yak (Bos grunniens; Hassanin & Douzery 1999; Verkaar et al. 2004) and the reciprocal monophyly of European aurochsen and taurine cattle (Troy et al. 2001; Edwards et al. 2007). Among taurine cattle, there was evidence for monophyly of haplogroups T1, T2 and T4 (table 2), but not for haplogroups T and T3.

View this table:
Table 2

Divergence time estimates for the bases of various bovine clades, made with and without assuming a proportion of invariant sites

Figure 2

Maximum clade credibility tree estimated using Bayesian analysis of 481 mitochondrial D-loop sequences from five bovine species, using the TIM+I+G substitution model. Numbers above nodes denote mean estimates of divergence times (kyr), with numbers below nodes showing the 95% HPD. The posterior probability of all the labelled nodes was 1.0. All branch lengths are scaled according to time; tips that do not reach the present represent ancient samples.

Estimates of the growth rate parameter for the five demographic models revealed significant evidence for population expansion in taurine and zebu cattle, and aurochsen, but not in bison or yaks. The result for aurochs is consistent with the findings of Edwards et al. (2007), where an exponential growth model provided a better fit to the data than a model assuming constant population size. Although our result for yak appears to differ from the qualitative result obtained by Finlay et al. (2007), we detect significant population expansions when the three major intraspecific clades are treated as separate populations (results not shown), with an average estimated age of 11 160 years. This treatment might be more appropriate due to the complex demographic history of yaks. There is evidence to suggest that Pleistocene climatic cycling between glacial and interglacial periods caused yaks to alternate between fragmented and contiguous interbreeding populations, resulting in large intra-population divergences within modern yaks (Guo et al. 2006). Therefore, the common ancestor of all yak lineages predates the most recent range expansion and does not correspond to the advent of domestication; instead, our mean age estimate for the common ancestor of yak mtDNA (approx. 75 kyr BP) coincides with the onset of the European Weichselian glaciation.

The coefficient of variation of rates among branches was 0.025 (95% HPD: 0.000–0.075) with invariant sites and 0.044 (95% HPD: 0.000–0.116) without invariant sites, suggesting that there is no significant departure from the assumption of a global molecular clock. When a proportion of invariant sites was assumed, the average substitution rate across the tree was estimated at 5.89×10−7 substitutions per site per year (95% HPD: 4.66×10−7–7.10×10−7 substitutions per site per year), equivalent to 58.9% per Myr. This dropped slightly to 5.68×10−7 substitutions per site per year (95% HPD: 4.31×10−7–7.02×10−7 substitutions per site per year) when the invariant sites' parameter was omitted.

Date estimates were made for the MRCAs of various bovine species and clades (table 2). The removal of the parameter for invariant sites had only a very small effect on most of the estimated dates. Collectively, they present a shorter time scale for bovine evolution than that suggested by previous studies. The estimated age of the root was 411 kyr (95% HPD: 269–564 kyr) with invariant sites and 374 kyr (95% HPD: 240–531 kyr) without invariant sites.

4. Discussion

Our results support a relatively contracted time frame for bovine evolution. For all major haplogroups of taurine and zebu cattle, basal divergence times were significantly more recent than 15 kyr ago, supporting the hypothesis of a domestication-associated population expansion. The molecular age estimates are consistent with the earliest evidence of domestic taurine cattle in the Near East, dating from the eighth millennium BC (Peters et al. 1999), and of zebu cattle in the fifth millennium BC (Bökönyi 1997). By contrast, the European aurochs seems to have experienced population growth concurrent with glacial retreat in Europe. The estimated age of the MRCA of all bison was very similar to that obtained by Shapiro et al. (2004), despite the differences in datasets and methodology.

The estimate for the age of the root of the tree is considerably younger than the 1 Myr BisonBos calibration point frequently used in the studies of cattle domestication (e.g. Bradley et al. 1996; Troy et al. 2001). The fossil record supports an older date for the BisonBos split, which falls between the first appearance of stem-lineage bovines at around 8.9 Myr (Barry et al. 2002) and the earliest representative of Bison dating from 2 to 2.5 Myr (Tedford et al. 1991). Overall, this suggests that it might not be appropriate to estimate dates for deeper, interspecific splits using aDNA-based calibrations at the tips of the tree, and that fossil-based calibrations at internal nodes might be preferable. Given the time dependency of rates (Ho et al. 2005), it might be advisable to employ a combination of the two calibration types to allow for higher rates towards the tips of the tree and lower rates in the basal branches.

Our internally calibrated phylogenetic analysis using multiple demographic models has yielded significant evidence of population expansion in a number of bovine species. The expansions in taurine and zebu cattle appear to coincide with domestication rather than post-glacial climatic changes. Further mitochondrial and nuclear data will serve to improve the precision of the current set of date estimates, allowing a more detailed insight into the evolutionary history of bovines.


This research was supported by the Australian Research Council, the Leverhulme Trust and the Royal Society. We are grateful to Alexei Drummond and Andrew Rambaut for making it possible to use multiple coalescent priors in their software Beast, Alan Cooper for assistance and three anonymous reviewers for their constructive feedback.


    • Received February 9, 2008.
    • Accepted March 28, 2008.


View Abstract