Diversification of Ramphastinae ( Aves , Ramphastidae ) prior to the Cretaceous / Tertiary boundary as shown by molecular clock of mtDNA sequences

Partial cytochrome b and 12S rDNA mitochondrial DNA sequences of eight representatives of the Ramphastidae family were analyzed. We applied the linearized tree method to identify sequences evolving at similar rates and estimated the divergence times among some of the taxa analyzed. After excluding Ramphastos tucanus and Capito dayi from our data set, the remaining taxa presented a constant rate of DNA substitution, and branch lengths could be re-estimated with a clock constraint using the maximum likelihood method. Branch lengths were calibrated assuming that Galliformes and Piciformes split around 100 million years ago (mya). Our results indicate that Ramphastinae, and probably Capitoninae, diverged from other Piciformes in the Late Cretaceous (~82 mya), suggesting that Piciformes is another avian order that survived the mass extinction event occurred 65 mya at the Cretaceous/Tertiary (K/T) boundary. The divergence times estimated within the Ramphastinae genera cover the period from the Middle Eocene (around 47 mya) through the Late Miocene (9.5 mya). Our estimate of divergence time is coincidental with the split of the African and the South American continents and other intense geologic activities and modifications of the areas which correspond to the current Neotropics. These events might have influenced the diversification of Ramphastinae in South America.


Introduction
The biodiversity of Neotropical regions ranks amongst the highest in the world.To better understand the mechanisms that generated such high levels of biodiversity, we should be able to identify when, where and how the separation between taxa occurred.
Whenever a representative fossil record from various strata of different ages is available, it makes it possible to trace back the origin of modern taxa and obtain reliable estimates of splitting events.However, fossil records are rather poor in the Neotropics as compared to those of other regions of the world.Although reliable dating of diversification events between taxa is available for plants through valuable information on pollen records and petrified wood, and even for some marine and aquatic fauna (e.g.Petri and Fulfaro, 1988), such information is very scarce for mammals and, except for some records on ancient birds, practically nonexistent for modern genera of birds.However, a few extinct species of Neotropical birds belonging to modern genera have been described for the Pleistocene (Alvarenga, 1993).
In the absence of fossil records, DNA sequencing and analysis is one of the available tools to estimate divergence times.This approach allows formulating hypotheses about the events related to taxa differentiation by comparing the estimated time elapsed since splitting of taxa with well-documented paleogeological and paleoclimatic events.For example, estimates for intergeneric diversification in the Neotropical avifauna are available for Psittacidae (Miyaki et al., 1998) and Cracidae (Pereira et al., 2002) and might provide valuable clues regarding the origin of taxa diversification in the Neotropics.
The family Ramphastidae occurs from Mexico through Argentina and comprises 55 species of Neotropical birds belonging to two subfamilies: Ramphastinae (toucans and aracaris, six genera) and Capitoninae (New World barbets, three genera) (Sibley, 1996).The interest in this group has been increasing exponentially, and several recent papers have reported on comparative studies on morphology, isozymes, hybridization and molecular sequences (e.g.Höfling, 1998;Castro et al., 2002;Barker and Lanyon, 2000).Barker and Lanyon (2000) analyzed 888 cytochrome b base-pairs from at least two species of each genus through a plethora of weighting schemes for maximum parsimony (MP) and proposed a phylogenetic hypothesis involving all Ramphastidae genera.Their results indicated Capitoninae (New and Old World barbets) as well as New World barbets to be paraphyletic, as suggested before (Prum, 1988;Sibley and Ahlquist, 1990), besides a very robust topology among the toucans and aracaris (Ramphastinae).Among these, the clade leading to Ramphastos was the first to split off, followed by Aulacorhynchus as a sister to ((Andigena, Selenidera), (Baillonius, Pteroglossus)).Maximum likelihood (ML) tests showed similar results, and most topological differences observed among ML and some of the MP weighting schemes might be related to stochastic errors due to small sequence sampling.The authors observed evolutionary rate variation among taxa, but they made no attempt to estimate the divergence time for any of the studied taxa.
Therefore, in the present study, we applied the linearized tree method (Takezaki et al., 1995) to a data set of partial mitochondrial cytochrome b (cyt b) and the small ribosomal subunit 12S rDNA sequences of all six valid Ramphastinae genera and one Capitoninae genus.We identified the taxa, which evolved at a rate that was significantly different from the others, and estimated the divergence times among those taxa, which evolved within the average rate.We also discuss a possible evolutionary scenario for their diversification, based on events of the history of Earth.

DNA sequencing
Total DNA was extracted from blood samples from seven Ramphastinae and one Capitoninae representatives (Table 1) as described elsewhere (Bruford et al., 1992).The DNA samples were used in PCR amplifications done in a final volume of 10.0 µL containing 20 ng of total DNA, polymerase buffer (Jeffreys et al., 1988), 1.0 µM of each primer, and 2.5 units of Taq DNA polymerase (Perkin Elmer).The primers employed to amplify and sequence a portion of the cyt b gene were L14841 and H15149 (Kocher et al., 1989), H16065 (Kornegay et al., 1993), reverse L (Cheng et al., 1994), L-internal, H-internal 2 and H-internal 3 (Miyaki et al., 1998).The 12S rDNA L and 12S rDNA H primers (Miyaki et al., 1998) were used to amplify a portion of the 12S rDNA gene.Thirty DNA amplification cycles were carried out, consisting of denaturation for 1 min at 95 °C, annealing for 30 s at 50 °C, and extension for 40 s at 72 °C each, using a TC2400 thermocycler (PE Applied Biosystems).An initial denaturation step of 5 min at 95 °C preceded the amplifications.A final extension step of 5 min at 72 °C was utilized.PCR products were separated on a 1% agarose gel and recovered in a 15% PEG8000 solution prepared with TBE buffer (Zhen and Swank, 1993).These   (Pereira et al., 2002).4: complete genome (Desjardins and Morais, 1990).
samples were then used in DNA sequencing reactions.A single PCR product of the expected size for cyt b and 12S rDNA genes was visualized on the agarose gel.Sequencing reactions were prepared using either the Dye Terminator Cycle Sequencing kit (PE Applied Biosystems) or the Thermo Sequenase Dye Terminator Cycle Sequencing kit (Amersham Life Science), according to the manufacturer's instructions.The samples were submitted to electrophoresis on an ABI PRISM 310 Genetic Analyzer (PE Applied Biosystems).Chromatograms of both strands of each gene were compared and the ambiguities were visually corrected using the Sequence Navigator software (version 1.0.1,PE Applied Biosystems; Parker, 1997).The GenBank accession numbers of all sequences obtained are listed in Table 1.This table also shows the accession numbers of Lanyon and Hall's (1994) Ramphastos tucanus sequence and Barker and Lanyon's (2000) Capito dayi sequence, which were included in our analyses.
The corresponding sequences from two other Piciformes (Melanerpes carolinus and Sphyrapicus varius) (Espinosa de los Monteros, 2000) and two Galliformes representatives (Crax blumenbachii and Gallus gallus) (Pereira et al., 2002;Desjardins and Morais, 1990) were included in the phylogenetic analysis, and their accession numbers are also shown in Table 1.The sequences were aligned in Clustal W (Thompson et al., 1994) and visually checked for misaligned positions in MacClade (version 3.08, Maddison and Maddison, 1992).

Phylogenetic analysis
Cyt b and 12S rDNA sequences were tested for multiple hits by plotting uncorrected distances against the number of transition and transversion for all pairwise comparisons.These analyses were carried out using MEGA (version 2.1, Kumar et al., 2001).All three codon positions of the cyt b sequences were included in the data set for the tree reconstruction.Sites containing missing information in both sequence sets and the gaps found in 12S rDNA sequences were excluded.Both cyt b and 12S rDNA sequences were concatenated in one final data set.
Tree reconstruction was performed by means of the branch and bound algorithm under maximum parsimony (MP) and maximum likelihood (ML) criteria, using PAUP* (version 4.0b8, Swofford, 2002) and selecting Crax blumenbachii and Gallus gallus as outgroups.To identify the best model of evolutionary changes for the ML analysis, the data set was subjected to the hierarchical likelihood ratio test (LRT, Huelsenbeck and Rannala, 1997), carried out using Modeltest (version 3.0, Posada and Crandall, 1998).
In this case, it tested whether the addition of parameters (such as unequal base frequency, rates of transition between purines and between pyrimidines, rates of transversions, proportion of invariable sites, and rate variation among sites) significantly improved the likelihood of the ML tree.Rate variation among sites was accounted for by three different assumptions.First, it was assumed that a proportion of sites is invariable (I) and that the variable sites evolve at the same rate; second, that all sites evolve according to a discrete gamma distribution (G); and finally that, differently from the first assumption, the variable sites are allowed to evolve according to G. The parameters obtained by Modeltest were used as user-defined parameters in PAUP, to speed up likelihood computation.Support for nodes in our reconstructions was obtained by jackknifing our data set in PAUP.Two sets of analysis were performed in this case, one assuming 30% and the other 50% deletion in each replicate, for 100 replicates.The purpose of the jackknife analysis was to measure the robustness of the Ramphastid phylogeny with a more reduced data set.
To determine whether the taxa included in our analysis evolved at the same rate, we performed a LRT.In this case, we tested whether likelihood values for trees with and without the assumption of a clock-like mode of evolution (Felsenstein, 1981) are significantly different.Unfortunately, the LRT was significant, suggesting that at least one sequence had an evolutionary rate that was significantly different from the others.
There are options available to determine divergence time in the absence of a general molecular clock.Some methods allow rates to differ among lineages (Sanderson, 1997;Thorne et al., 1998) and others accept rate variation at any point of the tree (Huelsenbeck et al., 2000).However, these studies present limitations and should be further investigated (Yoder and Yang, 2000).Alternatively, finding out which taxa evolve at a lower or higher rate than the average would help selecting and excluding the atypically evolving taxa.The analysis could then be performed with only those taxa, which evolve at similar rates (Takezaki et al., 1995;Rambaut and Bromham, 1998).The linearized tree method (Takezaki et al., 1995) consists of two different tests (two-cluster and branch-length), which seem to yield similar results and have been widely applied, for instance, in estimates of divergence time for the origin of orders of birds and mammals (Hedges et al., 1996), drosophilid species (Russo et al., 1995;Takezaki et al., 1995) and hominoids (Takezaki et al., 1995).Both tests require an outgroup in order to root the tree.The two-cluster test identifies statistical differences between the average substitution rates of two clusters joined by an interior node.This information allows eliminating those clusters that differ from the average root-to-tip distance.On the other hand, the branch-length test allows the identification of taxa evolving significantly faster or slower than the average.We decided to test for clock-like behavior through the branch-length tests implemented in the LinTree software (Takezaki et al., 1995).Here, we eliminated taxa evolving at significantly different rates and performed another round of the test, until no significantly different evolution rates Diversification of Ramphastinae prior to the Cretaceous/Tertiary were found among sequences (Takezaki et al., 1995), and molecular clock assumption was possible.
After excluding taxa with significantly different rates, we re-estimated branch lengths through maximum likelihood, assuming rate constancy among taxa and using the evolutionary model that best fit the data.We assumed that Galliformes diverged from other neognaths in the Cretaceous around 100 mya, as proposed by Paton et al. (2002).Standard errors for time estimates were obtained by non-parametric bootstrapping of our original data set.

Results
Our final alignment comprised 1160 sites (846 sites for cyt b, and 314 sites for 12S rDNA sequences).After the exclusion of the sites containing missing information and gaps, a total of 1128 sites were used for the phylogenetic reconstruction.No indels were present in the aligned cyt b sequences, and the translated amino acid sequences were highly similar to other vertebrate cyt b sequences available in the GenBank.The 12S rDNA sequences were also very similar to those described for other Piciformes and other birds deposited in the GenBank.Third codon positions of the cyt b sequences were saturated by multiple hits, as shown by comparisons between families, but not within families (data not shown).It has been shown that third codon positions can negatively affect parsimony analyses if saturated (Barker and Lanyon, 2000).However, these positions do contain phylogenetic information (Björklund, 1999;Edwards et al., 1991;Hastad and Björklund, 1998;Källersjö et al., 1999;Yoder and Yang, 2000) and were therefore maintained in our analyses.
Base composition among variable sites of the combined cyt b and 12S rDNA sequences was homogeneous among the taxa analyzed, except for Capito dayi (p = 0.001).The transition/transversion ratio estimated from the data set was equal to 2.21.Among the sites analyzed, 658 were invariable and 107 were parsimony informative.Mean uncorrected distances were 0.128 ± 0.031 substitutions per site among the Ramphastidae, 0.175 ± 0.016 substitutions per site among the Piciformes, and 0.186 ± 0.014 substitutions per site between Piciformes and Galliformes.
Tree topologies revealed by MP and ML analysis in the present work are presented in Figure 1 and are mostly congruent with those published by Barker and Lanyon (2000).The same tree topology was obtained when these genes were analyzed separately by either MP or ML.
LRTs revealed that the model which best fit the pattern of nucleotide substitutions found for these sequences was the transversional model (TVM, assuming unequal base frequency, different transition and transversion rates; Posada and Crandall, 1998), accounting for rate heterogeneity of substitutions among sites.Results are shown in Table 2.
We applied the LRT to test the molecular clock hypothesis and found that it did not support a clock-like behavior of nucleotide substitutions of the analyzed sequences.This result prevented the estimation of the divergence times among the Ramphastidae (log-likelihood without clock = -2966.38,log-likelihood with clock = -2979.43;χ 2 = 26.10,d. f. = 10, p = 0.0036).Thus, to select only the constantly evolving sequences, so that the molecular clock could be assumed, we performed the branch-length test.The results (Table 3) showed that R. tucanus and C. dayi have nucleotide substitution rates, which are lower and higher, respectively, than those of the remaining taxa (at α = 0.05, Figure 1).
Finally, after excluding R. tucanus and C. dayi from our data set, the remaining taxa presented similar rates of nucleotide substitutions, and thus branch lengths could be re-estimated with the clock constraint (Figure 2).Assuming that the Galliformes-Piciformes split occurred around 100 million years ago (mya, Paton et al., 2002), our estimates indicate that toucans and aracaris have diverged from other Piciformes around 82 mya (Table 4).Also, it is possible to estimate the generic diversification within the subfamily Ramphastinae occurred from the Middle Eocene (around 47 mya) to the Late Miocene (9.5 mya).The divergence time between Ramphastinae and Capitoninae (represented in this study by C. dayi) could not be estimated because the sequences of C. dayi did not evolve at the same rate as other Ramphastids.It should, however, be prior to the separation between toucans and aracaris, as Capitoninae is a sister group to Ramphastinae (Figure 1).

Discussion
In our study, partial cyt b and 12S rDNA sequences of eight Rhamphastidae representatives were aligned to correspondening sequences of two other Piciformes and two Galliformes.The evolutionary hypothesis provided in our study and that of Barker and Lanyon (2000) are very congruent showing the robustness of the classification at the generic level and therefore no further comments on the tree topologies will be made.However, we will focus on estimates of divergence times for the genera within the subfamily Rhamphastinae.
Estimation of divergence times based on molecular data is a common practice in phylogenetic studies (e.g.Arnaiz-Villena et al., 1999;García-Moreno and Silva, 1997;van Tuinen et al., 1998) and could offer insight into temporal mode (e.g.rates of diversification) or historical determinants of diversification, as geologically induced vicariance.Fossil data (e.g., Fleischer et al., 1998) and biogeographic events (e.g., García-Moreno and Mindell, 2000) can be used as calibration points for estimating molecular divergence times (Smith, 1998).However, for Ramphastidae, this information is very scarce (Feduccia, 1996;Sick, 1993).Three different approaches were available for the estimation of divergence times in our study: 1) the molecular estimate for the separation of Galloanserae from Neoaves 100 mya (Paton et al., 2002); 2) a mitochondrial DNA evolutionary rate of 2.0% per mya (Shields and Wilson, 1987;Kornegay et al., 1993;Klicka and Zink, 1997;Fleischer et al., 1998); and 3) the 0.5% rate for the transversions occurring at third codon positions per mya (Irwin et al., 1991).
We decided for the first option for the following reasons: (1) the molecular estimate is supported by a body of evidence as shown, for instance, by the fossil record (e.g., Stidham, 1998) Haddrath andBaker, 2001, Paton et al., 2002); (2) some groups seem to have evolved at rates other than the 2.0% sequence divergence per million years (Kessler and Avise, 1985); (3) using only transversions at the third codon positions would exclude 12S rDNA sequences and reduce our data set to fewer base pairs, since not all third codon positions of the cyt b sequences used here presented transversions, thus increasing the stochastic error in our estimates.Furthermore, the use of ancient calibration makes the estimation an interpolation rather than an extrapolation, which is inherently subject to greater errors.
The branch length tests indicated that the sequences of R. tucanus and C. dayi did not evolve at the same rate as the others.As pointed out by Barker and Lanyon (2000), there is non-stationarity in cyt b third positions of Capito.Incidentally, this base composition effect explains the apparent rate increase in Capito reported here.By excluding R. tucanus and C. dayi from our data set, we were unable to reject a clock-like behavior for the remaining taxa.
Divergence times, as estimated here, indicate that the subfamily Ramphastinae, and probably the Capitoninae, diverged from other Piciformes in the Late Cretaceous, well before the Cretaceous/Tertiary (K/T) mass extinction events of 65 mya.Indeed, Cooper and Penny (1997) estimated, based on fossil calibrations, that at least 22 lineages of modern birds were diversified prior to the K/T boundary.In their study, the authors did not include any member of Piciformes.Therefore, an important result of our study is the rise of the number of extant orders, which survived the K/T mass extinction.Although suggesting a cause for the diversification of Ramphastidae is a difficult task, some events of Earth history may help fulfill the endeavor.
The separation of Ramphastidae from other Piciformes, as estimated here, also coincides with the separation of Africa from South America.These landmasses began to fragment around 120-130 mya due to sea-floor spreading and, by 90 mya, they were completely separated (Salgado-Labouriau, 1994).Some evidence exists that the terrestrial vertebrates could probably move through these landmasses until 80-70 mya (Sibley and Ahlquist, 1990, and references therein).The range of divergence times within the subfamily Ramphastinae estimated here includes the period from the Middle Eocene to the Late Miocene.Those epochs are known to have had intense geologic activities that drastically changed the Earth climate and geography (Briggs, 1987;Salgado-Labouriau, 1994).One of the major changes includes the rising of the Andes (Potts and Behrensmeyer, 1992;Marshall and Sempere, 1993) that modified the climate and the river basins throughout the South American continent.A large sea transgression covering extensive areas of central South America has also been reported for that time (Petri and Fulfaro, 1988).These events might have influenced dispersal and diversification of living organisms, by changing or generating new habitats and even new ecosystems.For instance, the divergence time estimated here between the mountain-toucans, Andigena, and their lowland sister-taxon, Selenidera, was around 27 mya.This time corresponds to a major period of Andean uplift (Potts and Behrensmeyer, 1992;Marshall and Sempere, 1993), and we could assume that Andigena became adapted to the greater altitudes provided by the newly formed Cordillera of the Andes.
Most members of the subfamily Ramphastinae live in forests, but several Pteroglossus species live in drier environments (Haffer, 1974).These environments are currently represented by the cerrado (Brazilian savanna) and caatinga   (Brazilian semi-arid vegetation) biomes of Central and North-Eastern Brazil, respectively.We suggest that the drying out of central South America leading to the formation of the Brazilian savanna and the isolation of the Atlantic forest during the Miocene may have allowed the diversification of Pteroglossus from other toucans.
In conclusion, one of the significant results of our work was the inclusion of Ramphastidae among the lineages that survived the K/T extinction.Also, the uplift of the Andes during the Miocene and Eocene seems to have offered opportunities for the genera diversification of toucans, as well as of other Neotropical birds (Miyaki et al., 1998;Pereira et al., 2002) and of mammals (Cortés-Ortiz et al., 2003).These opportunities are related to geographic changes in South America, namely alterations in river basins, forest distribution, and formation of altitude and savanna-like habitats.A better understanding of how these processes influenced diversification within these groups should be provided by a more comprehensive study, increasing taxon sampling within Ramphastidae and other Neotropical groups (e.g.Cortés-Ortiz et al., 2003).

Figure 1 -
Figure 1 -Ramphastinae phylogeny estimated by cyt b and 12S rDNA sequences under the maximum likelihood criterion using the treebisection-reconnection algorithm on PAUP* 4.0b8.Numbers correspond to 30%-deletion/50%-deletion jackknife proportions.Ramphastidae members are underlined.The remaining taxa are Piciformes, except for Galliformes, Crax and Gallus, used as outgroup.Note the short and long branch lengths for R. tucanus and C. dayi, respectively.The scale bar corresponds to the expected proportion of nucleotide substitution per site.See results for likelihood value and parameters for this tree.
and also by mitochondrial and nuclear DNA sequences (e.g., van Tuinen and Hedges, 2001, Diversification of Ramphastinae prior to the Cretaceous/Tertiary 415

Figure 2 -
Figure 2 -Maximum likelihood clock tree under the transversional model with rate heterogeneity of DNA substitution among sites.Divergence times are indicated in a time scale below the phylogenetic hypothesis for those taxa evolving at a similar rate of DNA substitution.Gallus and Crax (Galliformes) were used as outgroup.Melanerpes and Sphyrapicus are members of Picidae, and the remaining taxa belong to Ramphastinae.

Table 1 -
Cyt b and 12S rDNA sequences of taxa analyzed in the present study.

Table 2 -
Hierarchical Likelihood Ratio Tests (HLRTs) performed for seven hypotheses of DNA substitution for our combined cyt b and 12S rDNA sequences.

Table 3 -
Branch length test under Tamura-Nei model of DNA substitution, assuming sequences to present rate heterogeneity of substitutions among sites (α = 0.18).

Table 4 -
(Paton et al., 2002)mates (in million years) for the Neotropical Ramphastinae and other Old World Piciformes, based on the divergence of Galliformes and other neognaths having occurred 100 million years ago(Paton et al., 2002).