Angiosperm diversification through time^†

Phylogenetic relationships

Phylogenetic trees used to estimate diversification rates were obtained from a data set including 265 genera belonging to 52 angiosperm orders (and unplaced isolated families), according to the Angiosperm Phylogeny Website (APW; 82). The APW is strongly based on the APG system (2; 3) and updates its contents with numerous recent publications. The orders sampled in this study represent over 94% of angiosperm living species. The gymnosperms Ginkgo biloba (ginkgophytes) and Gnetum gnemon (gnetophytes) were also included, and Ginkgo biloba was specified as the outgroup (Appendix S1, see Supplemental Data with the online version of this article). The data are the concatenated nucleotide sequences of three plastid genes (atpB, rbcL, and matK) and two nuclear markers (18S nuclear ribosomal DNA [nrDNA] and 26S nrDNA). Sequences were obtained through a bioinformatic search in GenBank, which aimed to balance a comprehensive representation of angiosperm orders with the largest possible number of molecular markers for all included taxa. The matrix was complete in that all genes were available for all taxa, although in several cases, a single genus was represented by different species. Sequence alignment for each marker was achieved with the program MUSCLE (21), followed by manual refinements with the program BioEdit v5.0.6 (34). The sequences of each marker were subsequently concatenated in a single data set (http://treebase.org, TreeBASE accession SN4019).

An examination of parameter values of best fitting models for individual molecular markers and codon position partitions (for atpB, rbcL and matK), obtained with the Akaike information criterion (AIC) implemented in the program Modeltest (61; 60), indicated that the data could be appropriately divided into four partitions: (1) first and second positions of atpB and rbcL, (2) third positions of atpB and rbcL, (3) matK, and (4) 18S and 26S. Phylogenetic relationships among the 265 angiosperm genera were estimated with Bayesian analysis using MrBayes v3.1.2p (36). Independent models with unlinked parameters were applied to the four data partitions, implementing variable rate priors. Two independent Metropolis coupled-Markov chain Monte Carlo (MC³) runs of 5 × 10⁶ generations, each consisting of four incrementally heated chains (temp = 0.2), were conducted, sampling one tree every 200 generations. After examination of generation vs. likelihood plots, trees corresponding to the initial 500000 generations (10%) were discarded as burn-in. Post burn-in sampled topologies were summarized as a 50% majority rule tree to obtain the posterior probability (PP) credibility interval of each clade.

Age estimation

Ages of clades were estimated with penalized likelihood (PL; 69, 15), a molecular-based semiparametric method that incorporates among-lineage rate heterogeneity and can use fossil information as auxiliary in divergence time estimation. We relied on two considerations to select a phylogenetic tree to estimate dates. First, avoid an incompletely resolved tree (e.g., the 50% majority rule consensus), to circumvent possible complications in the cross validation procedure (described later), probably the most computationally difficult step in PL dating (71, pp. 15–16). Second, use a phylogram (i.e., topology with branch lengths) in which branch lengths were obtained through the optimization of best-fitting model parameters for each data partition. Hence, we selected the topology with highest PP found in the two Markov chains as a phylogenetic working hypothesis and used all the phylograms that have a topology identical to it to estimate ages. The topology with highest PP was identified with MrBayes, and topologically identical phylograms were found and extracted with the program PAUP* version 4.0b10 (84). Because of technical conventions for optimizing branch lengths around the root node of a tree (e.g., 71, pp. 23–24), it became necessary to remove Ginkgo biloba, the outgroup used during phylogeny estimation, from dating analyses. The divergence between Gnetum gnemon and angiosperms became the new root node.

Phylograms were temporally calibrated by fixing the age of the root node, which represents the crown group of seed plants, at 350 Myr. This age was selected because it is younger than the Famennian (Upper Devonian; 67) age of the oldest known fossil seeds (Elkinsia polymorpha), and older that the Namurian (Lower-Upper Carboniferous; 85) age of the oldest presumed crown group seed plants (Cordaitales). It also corresponds approximately to the mean age for the crown group of seed plants obtained in different relaxed molecular clock analyses for a representation of vascular plant lineages (S. Magallón, unpublished results). Forty-nine nodes within angiosperms were constrained with minimal ages (minage) derived from a critical examination of the angiosperm fossil record (Fig. 1, Appendix S2, see Supplemental Data with the online version of this article).

In addition to the minage constraints, two maximal age (maxage) constraints were alternatively applied. In the first case, here referred to as “relaxed” dating, a maxage of 125 Myr was applied to the eudicot crown node. The appearance of tricolpate pollen grains in late Barremian-early Aptian (ca. 125 Myr) sediments (39; 17; 37; 31), has been regarded among the best indications in the fossil record of the origin of a biological lineage. There is an exact correspondence between the presence of a distinctive morphological attribute, which was abundantly produced and became easily fossilized (i.e., tricolpate pollen), and membership to a monophyletic group (i.e., eudicots, or tricolpate angiosperms). The age of the eudicot crown node appears to be constrained narrowly around 125 Myr, given the Barremian-Aptian age of the oldest tricolpate pollen, presumably produced by eudicot stem lineage representatives, and the Barremian-Aptian report of putative eudicot grown group members (47; 31). Previous studies (e.g., 79; 1) have used the earliest record of tricolpate pollen as a definitive temporal landmark for the origin of eudicots.

In the second case, referred to as “constrained” dating, a maxage of 130 Myr was applied to the angiosperm crown node. This maxage is based on the oldest report of fossil angiosperm pollen, typically characterized by a thin and granular endexine (31), from Valanginian-Hauterivian (ca. 136.4 Myr) sediments (38; 40; 37; 9). Assuming that these early pollen grains were produced by angiosperm stem lineage representatives, we speculate that the origin of the angiosperm crown group may have occurred shortly afterwards, by the Hauterivian-Barremian boundary (130 Myr). The age of the angiosperm crown group is limited by a lower (younger) bound of approximately 125 Myr, corresponding to the late Barremian-early Aptian age of unequivocal angiosperm crown group members (Chloranthaceae, 31; Winteraceae, 20; 18; eudicots, described earlier). The stratigraphic position of fossils used for tree calibration and, as minage and maxage constraints, was transformed into absolute ages (Myr) using the upper (younger) bound of the interval, based on the stratigraphic time scale of 33.

The Langley–Fitch (46) test of the molecular clock, as implemented in the program r8s version 1.71 (70, 15) was conducted. Penalized likelihood, implemented in r8s, requires a user-defined parameter to specify the level of molecular rate smoothing to be implemented in dating analysis. To identify the smoothing magnitude (λ) that best describes the available data, we used a cross validation procedure that calculates the predictive error associated to molecular rate and time estimates across the full tree, derived from sequentially removing minage and maxage constraints (71). Cross validations were performed for relaxed and constrained maxage implementations, using one randomly selected phylogram among those with highest PP topology. Each cross validation tested 16 smoothing magnitudes ranging from logλ₁₀ = −2 to 5.5 at 0.5 intervals, which comprise a broad spectrum of substitution regimes. Penalized likelihood dating was conducted on all the phylograms with highest PP topology, using relaxed and constrained maxage constraints and implementing in each case the optimal smoothing magnitude found in the respective cross validation. Penalized likelihood analyses used a TN algorithm with bound constraints with five initial starts and three perturbed restarts with perturbations of magnitude 0.05 in random directions. The mean age, standard deviation, and 95% confidence interval of every node in the tree were derived from the point estimates of age in each of the phylograms with highest PP topology.

Diversification rates

Diversification rates were calculated using method-of-moments estimators (66) in the context of a birth-and-death model (44) that consider the species diversity and age of a clade. These estimators provide absolute estimates of the rate of diversification of a clade, they are conditional on the survival of the clade to a given time t, in this case, the present, and they can differentially estimate the diversification rate of a stem clade or of a crown clade (19; 53). These conditional estimators of absolute diversification rates correspond to eqs. 6 and 7 in 53: for stem groups and for crown groups. The relative extinction rate (ε) is defined as ε = μ/λ. The variable t corresponds to a time after the origin of the clade, here the present, and n is the standing species diversity of the clade at time t. Because absolute speciation and extinction rates for angiosperms and angiosperm clades are unknown, diversification rates were estimated assuming that the relative extinction rate is bounded within ε = 0.0, which implies no extinction, and ε = 0.9, which implies a very high relative extinction rate. Whereas ε = 0.0 represents an absolute lower bound for the relative extinction rate, the selection of ε = 0.9 as an upper bound is arbitrary. It is nevertheless justified by the observation that as ε approaches 1, the magnitudes of the rates of speciation (λ) and extinction (μ) increase rapidly, exceeding values of one event per million year. These values approximately represent the empirical maximum estimated from real data for a variety of animal taxa (80; 41). Also, for large clades, values of ε larger than 0.9 correspond to clades having less than a 10% chance of surviving to the present, which in the case of angiosperms, with a living diversity of about 270000 species, seems unlikely (53). Therefore, we estimated the diversification rate of each clade considering ε = 0.0 and 0.9. Standing species diversity in angiosperm orders was obtained from the APW (82).

In addition to angiosperm orders, absolute diversification rates were estimated for nonnested crown clades appearing every 10 Myr since the onset of angiosperm crown group diversification. This approach aimed to evaluate diversification rate from a standpoint that relies on the temporal distribution of branching events on a tree. The identification of nonnested crown nodes was based on the timing of angiosperm lineage diversification according to constrained dating. To correctly quantify the number of living species derived from a node in the tree, we needed to consider the species diversity of orders (and unplaced families) that were not sampled in the phylogenetic analysis. Unsampled orders were intercalated in the dated tree according to their position in the order-level phylogeny in the APW (82; Fig. 1). The stem group age of every intercalated order was estimated as the midpoint between the ages of dated nodes immediately above (younger) and below (older) it. Orders (or families) of unresolved phylogenetic position in the APW (82) were placed in a polytomy (Fig. 1).

We estimated diversification rates of angiosperm orders using relaxed and constrained dates; their stem group age (for 72 orders: 52 sampled, 20 inserted) and their crown group age (for 41 orders), with a relative extinction rate (ε) of 0.0 and of 0.9. Diversification rates of nonnested crown clades were also estimated using ε = 0.0 and 0.9. Estimated diversification rates were plotted against the age of the respective clade (Figs. 2, 3).

We attempted to implement topology-based, whole-tree statistical tests for detecting significantly different diversification rates in the angiosperm tree and to identify the branches of the angiosperm tree in which diversification rate shifts have occurred (Symmetree program; 11). Nevertheless, the tests could not be conducted due to the large number of species encompassed in several terminals of the angiosperm order-level tree (B. R. Moore, U.C. Berkeley, personal communication).

Phylogenetic relationships

The concatenated sequences of the five molecular markers (atpB, rbcL, matK, 18S nrDNA, and 26S nrDNA) are 9789 base pairs (bp) in length and are available for all the genera in the data set. The best-fit model for each of the four data partitions (i.e., 1^st and 2^nd codon positions of rbcL and atpB; 3^rd positions of rbcL and atpB; matK; and 18S-26S nrDNA) included six nucleotide substitution parameters, among-site rate variation, and a proportion of invariable sites (GTR+I+Γ). After 5 × 10⁶ generations, the two independent Markov chains were estimated to have converged (standard split frequencies < 0.01). Plots of generation number vs. likelihood value indicated that in both runs, likelihood values stabilized approximately after 200000 generations; however, the trees sampled in the initial 500000 generations (10% of the total) were excluded as burn-in. The 50% MR consensus of post burn-in trees (available in TreeBASE SN4019) contains five unresolved polytomies: within Alismatales, among major core eudicot lineages, within Saxifragales, and within Ericales (two polytomies).

The topology with the highest PP was selected as a working phylogenetic hypothesis for molecular dating. An order-level summary is shown in Fig. 1, and the full tree is available in TreeBASE (SN4019). Most relationships are in agreement with independent assessments (e.g., 78; 82), and relatively unusual relationships, for example, the sister group relationship between Amborellales and Nymphaeales, are weakly supported. Unless otherwise indicated, the following relationships are supported by ≥0.95 PP. The deepest split within angiosperms separates a branch including Amborellales and Nymphaeales (0.74 PP) from all other members of the clade. Austrobaileyales is the sister to core angiosperms, which are divided into a branch that includes Chloranthales plus monocots (0.57 PP) and eumagnoliids, and a branch that includes Ceratophyllales and eudicots. Ranunculales is the sister to all other eudicots, followed by a branch that includes Sabiales and Proteales (0.73 PP). Buxales is the sister to core eudicots. The 50% majority consensus tree contains a trichotomy involving Gunnerales, Dilleniales, and a clade that includes all other core eudicot lineages. In the highest PP topology, Gunnerales and Dilleniales are sister taxa (<50% PP). The deepest split within the clade that includes all other core eudicot lineages separates a clade (0.86 PP) that includes Saxifragales, Vitales, and rosids, and another clade that includes Santalales, Berberidopsidales, Caryophyllales, and asterids (Fig. 1).

Age estimation

The age and phylogenetic position of the 49 minage constraints are shown in Fig. 1 and described in Appendix S2 (see Supplemental Data with the online version of this article). Among the phylograms sampled by the two Markov chains, 46 are topologically identical to the topology with highest PP. Dating was conducted on these 46 phylograms.

The Langley–Fitch test of the molecular clock indicated that relaxed and constrained phylograms departed significantly from rate constancy (P = 0 in all tests). Not surprisingly given the number of taxa and the branch length heterogeneity in the phylogram, the relaxed and constrained cross validation procedures returned several failed optimizations. Nevertheless, both showed relatively small differences in the raw error associated to the range of tested smoothing magnitudes (λ). A general positive relation between raw error and smoothing magnitude was found, except for the highest smoothing magnitudes in the relaxed cross validation (results available from the authors). In both cross validations, the lowest smoothing magnitude that was tested (which passed the optimization) has the smallest associated raw error, and hence, relaxed and constrained dating analyses were conducted implementing a smoothing parameter of log λ₁₀ = −2 (λ = 0.01).

Ages of nodes derived from relaxed dating were substantially older than those obtained in the constrained analysis, except within eudicots, where ages are similar in both analyses. Both sets of ages were used to estimate diversification rates. The mean age and 95% confidence interval estimated for the stem group and crown group of angiosperm orders according to relaxed and constrained dating are shown in Table 1. The chronogram resulting from constrained dating is shown in Fig. 1. The relaxed chronogram is available from the authors.

Table 1. Species diversity, age, and diversification rate at two relative extinction rates (ε) for angiosperm orders. The two subheadings above each column describe the first and second row, respectively, for each order listing.

Nonnested crown group clades

Twenty orders (and unplaced families) were inserted in the constrained dated tree based on their position in the APW order-level tree (82) and a stem group age intermediate between that of immediately deeper and shallower dated nodes (Fig. 1). In many cases, the bounding nodes are separated by a short time gap; thus, the possibility of a large error in the age assigned to the intercalated nodes is small. Diversification rates were estimated for the inserted orders (Table 1); however, these should be regarded as especially tentative due to the unavailability of a direct estimate of age, and in some cases, uncertainty in phylogenetic position.

Thirteen nonnested crown clades appearing every 10 Myr and in six intervals from 130 to 70.01 Myr were identified (including angiosperms as a whole; Table 2). The three intervals between 120 and 90.01 Myr gave rise to the largest number of nonnested crown clades (three, two, and five, respectively), while the first (130–120.01) and two last intervals (90–70.01) each gave rise to one (Fig. 1, Table 2). The phylogenetic equivalence and species diversity (including inserted orders) of nonnested crown clades are shown in Table 2.

Diversification rates

Diversification rates calculated for angiosperm orders based on their stem group (72 orders: 52 sampled and 20 inserted) and crown group age (41 orders), derived from ages obtained with relaxed and constrained dating, and relative extinction rate (ε) of 0.0 and 0.9 are shown in Table 1. The absolute (crown group) diversification rates estimated for angiosperms as a whole are 0.0420 Myr⁻¹ with ε = 0.9, and 0.0489 Myr⁻¹ with ε = 0.0 using relaxed dating, and 0.0781 Myr⁻¹ with ε = 0.9, and 0.0909 Myr⁻¹ with ε = 0.0 using constrained dating. Angiosperm orders with the highest stem group diversification rates are Lamiales (0.0968–0.1257 Myr⁻¹), Gentianales (0.0928–0.1214 Myr⁻¹), Asterales (0.0834–0.1081 Myr⁻¹) and Solanales (0.0777–0.1074 Myr⁻¹). Orders with the lowest stem group rates are Amborellales (0 Myr⁻¹), Acorales (0.0015–0.0109 Myr⁻¹), Ceratophyllales (0.0020–0.014 Myr⁻¹), and Berberidopsidales (0.0023–0.0123 Myr⁻¹). Among the orders for which crown group age was available, Malvales was identified as having the highest diversification rate (0.1877–0.2366 Myr⁻¹), followed somewhat distantly by Lamiales (0.1225–0.1490 Myr⁻¹), Brassicales (0.0927–0.1182 Myr⁻¹), Asterales (0.0928–0.1127 Myr⁻¹), and Sapindales (0.0879–0.1117 Myr⁻¹). Crown group clades with the lowest rates are Berberidopsidales (0.0023–0.0078 Myr⁻¹), Austrobaileyales (0.0115–0.0313 Myr⁻¹), Canellales (0.0134–0.0315 Myr⁻¹), and Chloranthales (0.0135–0.0290 Myr⁻¹).

All else being equal, diversification rates derived from constrained dating were typically higher than those derived from relaxed dating, and those obtained with ε = 0.0 were higher than those obtained with ε = 0.9. In clades for which stem group and crown group diversification rate were estimated, the two rates were of similar magnitude. A general positive relationship between the stem group and crown group diversification rate for any given clade was found, and a tendency for one rate to be larger than the other was not observed.

Diversification rates of the 13 nonnested crown clades appearing every 10 Myr implementing ε = 0.0 and ε = 0.9, are shown in Table 2. The diversification rates estimated for nonnested crown clades appearing between 130 and 90.01 Myr are similar to those of angiosperms as a whole. Among these, only the clade corresponding to the most recent common ancestor (MRCA) of Pandanales and Dioscoreales (120–110.01 Myr interval) and the MRCA of Commelinales and Zingiberales (100–90.01 Myr interval), have somewhat lower rates (0.0457–0.0596 Myr⁻¹ and 0.0568–0.0734 Myr⁻¹, respectively; Table 2, Fig. 3). However, the diversification rates of the two youngest nonnested crown clades are substantially higher than all the rest (Table 2, Fig. 3). Diversification rates for the MRCA of Vahliaceae and Solanales (90–80.01 Myr interval) range from 0.1008 to 0.1207 Myr⁻¹, and for the MRCA of Lamiales and Solanales (80–70–01 Myr interval) from 0.0995 to 0.1202 Myr⁻¹. As observed for angiosperm orders, diversification rates obtained with ε = 0.0 were higher than those obtained with ε = 0.9.

Diversification through time

The relationships of diversification rate and clade age show that the highest diversification rates are found among younger orders, whether considering relaxed or constrained dating, stem or crown group diversification rates, or high or low relative extinction rates (Fig. 2). The diversification rate of nonnested crown clades is approximately constant, except for the two youngest nonnested crown clades, which exhibit appreciably higher diversification rates (Fig. 3).

Penalized likelihood dating

Penalized likelihood analyses included a comprehensive representation of angiosperms at the order level and incorporated a substantial amount of fossil-derived information. The 130 Myr maxage constraint to the angiosperm node in constrained dating exerted a determinant influence on age estimates across the tree. Constrained dates were much younger than relaxed dates; however, in some cases, they were substantially older than the earliest fossil record of particular angiosperm lineages (Fig. 1). Other studies (e.g., 87; 7; 54) have obtained ages for angiosperms or angiosperm clades that, as in the relaxed analysis, are substantially older than their oldest fossils. Here, we followed the approach of 74, who performed alternative dating analyses by fixing or unfixing the age of angiosperms at 132 Myr. These authors also found that ages obtained in the “unfixed” analysis were much older than those in the “fixed” analysis and the angiosperm fossil record.

In both relaxed and constrained analyses, the timing of appearance of orders is more or less continuously distributed through a delimited time interval. In relaxed dating, this interval spans approximately from 236 Myr (Middle Triassic) to 77 Myr (mid-Campanian, Lower Cretaceous) for stem groups, and from 203 Myr (Upper Triassic) to 34 Myr (ca. Eocene-Oligocene boundary) for crown groups; Fig. 2A, B). In constrained dating, the appearance of order stem groups spans approximately from 130 Myr (Hauterivian-Barremian boundary) to 77 Myr (mid-Campanian), and of order crown groups from 126 Myr (ca. Barremian-Aptian boundary) to 34 Myr (ca. Eocene-Oligocene boundary; Figs. 1, 2C, D).

Diversification rates

Diversification rates estimated here for angiosperms as a whole, using constrained dates, are slightly higher than those obtained previously with the same estimators (0.0767 Myr⁻¹ with ε = 0.9, and 0.0893 Myr⁻¹ with ε = 0.0; 53), as a result of the slightly older age for the angiosperm crown node (132 Myr vs. 130 Myr in constrained dating) and the slightly smaller species diversity (262196 spp. vs. 269323 spp. here) used in the previous study. Diversification rates of angiosperm orders vary extensively, for example, from zero to 0.13 Myr⁻¹ for stem groups, and from 0.002 to 0.15 Myr⁻¹ (or up to 0.24 Myr⁻¹ for Malvales, see ahead) for crown groups. Angiosperm orders with the highest relaxed and constrained stem group and crown group diversification rates belong to the lamid clade (Lamiales, Gentianales, Solanales, and Boraginaceae), the campanulid clade (Asterales, Apiales, and Dipsacales), and Ericales, within the asterids; to the malvid (Malvales, Brassicales, and Sapindales) and fabid (Malpighiales, Fabales, and Rosales) clades, within the rosids, and to the commelinid monocots (Poales). All these clades are nested highly within angiosperms. Notably, the rate estimated for Asparagales, which contains the megadiverse family Orchidaceae (orchids), is not among the very highest among angiosperms, possibly as a consequence of its very old stem and crown ages.

Diversification rate estimation hinges on the species richness and the age of a clade. Given constant species richness, an incorrect clade age will result in underestimated or overestimated diversification rates, if the used age is older or younger, respectively, than the true age. Among the many factors that influence molecular clock dating (e.g., 68), one that is crucial for accurately estimating a clade's crown age is the representation of the deepest phylogenetic split within the clade. If this bifurcation is not represented, but instead, only members of one of the branches stemming from it are included, the estimated age will be younger than the crown group age. Thus, just as fossil first appearances, albeit for different reasons, molecular estimates of crown group ages are minimal estimates, unless the deepest split within the clade is sampled.

Considering the small taxon sampling within angiosperm orders used in this study, it is very likely that the deepest split of each angiosperm order was not sampled. Hence, the “crown group ages” presented here may correspond to divergences within orders, and the order's true crown group age is older, maybe substantially so. Diversification rates derived from such artificially young ages are overestimates of the order's true crown group diversification rate. Perhaps, this artifact explains some of the unexpectedly high crown group diversification rates obtained for some angiosperm orders, namely, the Malvales (See Results, Table 1). With a standing diversity of nearly 6100 species and an estimated crown group age of 33.9 Myr, the crown group diversification rates of Malvales (0.1877 and 0.2366 Myr⁻¹ for ε = 0.9 and ε = 0.0, respectively) far exceed those estimated for any other angiosperm order. Nonetheless, although we suspect that crown group diversification rates of Malvales (and possibly also Brassicales and Sapindales) are inflated to some extent as a result of artificially young crown group ages, these clades are probably among those with highest diversification rates among angiosperms, as shown by their high stem group diversification rates (Table 1).

Diversification rates estimated for orders that were a posteriori inserted in the phylogeny are especially tentative because their stem group age was arbitrarily considered as the midpoint between the age of the immediately younger and older dated bounding nodes (Fig. 1). If the true age of an inserted order is older than the assigned age, the estimated diversification rate will be an overestimate of the true rate, and vice versa. Fortunately, for all inserted orders, the difference in age of the two bounding nodes is small, and thus, the range of possible ages that correspond to each inserted order is narrow. We expect that the implicit error in the age of inserted orders has a small effect on the estimated rates.

The higher diversification rates resulting from ages derived from constrained analyses, compared to those from relaxed analyses for any given node, simply reflect that, all else being equal, less time was available to produce the same number of species given constrained ages. The higher diversification rates obtained with ε = 0.0 with respect to ε = 0.9 may appear somewhat counterintuitive, but can be explained with the following reasoning: When ε = 0.0, the extinction rate (μ) is zero, and the rate of diversification (r = λ – μ) is equal to the rate of speciation (λ). As ε = μ/λ approaches one, the difference between the magnitude of λ and the magnitude of μ becomes smaller, and therefore, r decreases in magnitude.

Diversification rates derived from stem group age and from crown group age were similar in most of the 41 orders for which both could be estimated. Such correlation should not be necessarily expected. For example, a substantial difference between stem and crown diversification rates would occur in orders with an ancient stem group age (time of lineage differentiation) and a young crown group age (time of diversification to extant species). This study also revealed a general positive correlation between species diversity and diversification rate (results not shown). Although this correlation does not seem surprising, it is not obligatory. Species diversity and diversification rate may be decoupled, for example, in ancient lineages with a low diversification rate, which have accumulated a large number of species through time, or in young lineages with a high diversification rate, which have not had sufficient time to accumulate a high species diversity.

Delimitation of orders and the use of nonnested crown clades

The rate of diversification estimated for an order depends on its age and its species diversity, which in turn are contingent on the delimitation of the order. The taxonomic decisions used to define particular clades as “orders” are confounding factors on rate estimation; however, the effects of different taxonomic delimitations on estimated rates are not straightforward. Suppose that an alternative taxonomic system lumps into a single order two or more closely related APW-defined orders. The new order will have a greater species diversity, causing diversification rate to increase, but it will also have an older age, causing diversification rate to decrease. The extent to which different taxonomic delimitations cause substantial differences in diversification rate, this is, the extent to which change in species diversity is compensated by change in age, depends on the distribution of species diversity among the lumped or splitted subclades, their relative phylogenetic position, and their relative age. The interactions among these three factors, as well as the magnitude of the effect of alternative taxonomic delimitations, require explicit examination.

The use of nonnested clades that appear at regular time intervals represents an alternative to explore angiosperm diversification that considers branching events on the tree. The number and composition of nonnested crown clades identified here was greatly determined by the arbitrary division of time into 10 Myr intervals, and, significantly, by the number of terminal taxa included in the tree. Ideally, identification of nonnested crown clades through time should be based on a complete, dated species-level phylogeny. In spite of these important caveats, the use of clades delimited by branching through time provided an assessment of the long-term trends of diversification among the deepest angiosperm branches that is independent from taxonomic delimitations. Most of the identified nonnested crown clades appeared between 120 and 95 Myr ago (approximately mid-Aptian to mid-Cenomanian; Table 2, Fig. 3). Only two appeared later, between 85 and 75 Myr ago (approximately Santonian to mid-Maastrichtian; Table 2, Fig. 3). The rates of diversification of the two youngest nonnested crown clades are appreciably higher than those of the rest.

Increase of angiosperm diversification rates through time

Diversification through time plots show that the highest rates of diversification are found among the younger angiosperm orders. Nevertheless, age is not an effective predictor of relative diversification rate because some young orders have a moderate or low diversification rates (Fig. 2). Diversification through time plots for nonnested crown clades show approximately constant rates through angiosperm history, except for the two youngest nonnested crown clades, which have substantially higher rates (Fig. 3).

It could be argued that the increase in the rate of diversification rate through time represents the initial species accumulation of young clades during an exponential diversification process. However, the estimators used here are explicitly modeled as a stochastic birth-and-death process (4) where the speciation rate and the extinction rate are constant through time and lead to exponentially increasing or decreasing diversity (53). The diversification rate of a clade going through the initial phase of explosive species accumulation in an exponential process should not be detected as a higher net diversification rate because the estimators are expected to account for the initial explosive phase.

An equally serious concern regards the possibility that the highest diversification rates found among the youngest clades is an artifact resulting from retrospective inferences about diversification rates based exclusively on living diversity (63; 58). The accumulation of species diversity through time can be inferred from a molecular phylogeny by extracting branching times (e.g., with molecular clock methods) and plotting the cumulative increase in the number of species through the history of the clade (i.e., a lineage-through-time plot). 58 described the expected and observed trends of lineage-through-time plots for a clade diversifying under an exponential birth-and-death model in which speciation is larger than extinction, and both remain constant through time. The lineage-through-time plot derived from all living and extinct species belonging to a clade is a straight line with a positive slope. But, if only living species are considered, the lineage-through-time plot will have a terminal upward turn, caused by the presence of recently originated species that have not yet had time to go extinct. This terminal upward turn, termed the “pull of the present” (58), corresponds only to the speciation component of the diversification process and becomes more pronounced as the magnitude of extinction becomes greater (58). An explanation of the higher diversification rates exhibited by younger clades as an effect of the “pull of the present” should entail a lineage-through-time plot (and not the diversification-through-time plots presented here) based on species-level dated phylogenies. Furthermore, our finding that younger clades have higher diversification rates cannot be explained as “pull of the present” because the estimators of diversification rate used here account for the possibility of extinct, unobserved diversity, by being conditional on the survival of the clade to a given time after its origin (53), in this case, the present. The “pull of the present” (58) should be distinguished from the related but different idea of the “pull of the Recent” (62), which explains the apparent increase in fossil diversity as the Recent is approached as the result of built-in factors of the fossil record and the conventional methods of analyzing it (62). The main cause of the “pull of the Recent” is the extension of stratigraphic ranges of fossil taxa by the more compete sampling of the Recent biota (62; 42).

Diversification rates among ancestral angiosperms

The oldest living angiosperm lineages have some of the lowest detected diversification rates. Amborellales, with a single living species, has a net diversification rate of zero because the number of extant species (one) does not represent an increase from the number of species present in the lineage immediately after its phylogenetic differentiation. Diversification rates estimated for Nymphaeales, Austrobaileyales, and Chloranthales are also among the lowest among angiosperms, an inescapable consequence of their small standing diversity and old age.

One outstanding question is whether ancestral angiosperm lineages have contained few species through their history or whether their living species are merely a remnant of greater historical diversity. Diversification estimators used in this study are based on a birth-and-death process regulated by rates of speciation and extinction that are constant through time (i.e., time-homogeneous). The extent to which a time-homogeneous model adequately describes the process of diversification through time of different angiosperm lineages requires independent evaluation and is potentially a major limitation of this study. The existence of ancient and species-poor lineages is an improbable outcome under a variety of diversification conditions (83). Explaining them under time-homogeneous models usually requires that speciation and extinction rates were both very low and that a large number of independent lineages underwent the process simultaneously (83). An alternative explanation involves a time-heterogeneous model in which speciation was much higher than extinction during a short period after the onset of the diversification process, followed by a long period during which speciation and extinction were approximately equal (83).

If the diversification of ancestral living angiosperms corresponds more closely to a time-homogeneous process, it is likely that they have included few species through their history and that their extant diversity is the outcome of slow species accumulation through time. If their diversification corresponds more closely to a time-heterogeneous process, then ancestral angiosperm lineages are depauperate survivors of a much larger, currently extinct diversity. The possible representation of Amborellales (29, 31), Nymphaeales (30), Austrobaileyales (32; 31), and Chloranthales (28, 31; 23) in the Early Cretaceous, among the earliest angiosperm whose phylogenetic affinity can be securely identified (31), suggests a more extensive past diversity. If living ancestral angiosperm lineages encompassed a substantially larger species diversity in the past, their diversification may have included an explosive speciation phase at the onset of their evolution, which may later have been limited by competing lineages (most likely other angiosperms, e.g., 49; 51) or by resource availability (75).

The living ancestral angiosperm lineages represent a negligible proportion of extant angiosperm diversity. The nature of their evolutionary diversification lies somewhere between long-lasting perseverance under low speciation and extinction, and (more likely) survival from a much greater past diversity (discussed earlier). These lineages, however, represent crucially important elements of extant angiosperm composition because of their unique potential to document the vegetative and reproductive traits, physiological capabilities, and ecological role of angiosperms during the initial phases of their evolution.

Concluding remarks

This study provides estimates of the diversification rate of angiosperms as a whole and of a set of orders that represent a substantial proportion of living angiosperm species. Ages used to estimate diversification rates were obtained through molecular dating analyses that accounted for molecular rate heterogeneity, incorporated fossil information about minimal ages of lineages, and included or excluded a maximal constraint to angiosperm age. Diversification rate estimators are based on a time-homogeneous, birth-and-death exponential model and are contingent on lineage survival to the present.

The temporal distribution of the stem group and crown group diversification rates of angiosperm orders document substantial changes through angiosperm history, with a trend toward increasing rates among younger orders. However, orders with low rates occur throughout angiosperm history. The youngest nonnested crown clades were also found to have the highest rates. The high diversification rates found among the younger orders or crown clades cannot be attributed to the explosive phase of species accumulation during exponential growth, nor to a “pull of the present” effect, because the used estimators are based on an exponential diversification process, and they are contingent on the survival of the lineage to the present.

The extraordinary species richness of angiosperms is found to have mixed origins. Slightly less than half (44.5%) of living species belong to orders with low to moderate diversification rates that appeared, according to constrained dating, between 130 and 102 Mya, approximately from the Barremian to the uppermost Albian, in the Lower Cretaceous. Slightly over half of living species (55.5%), however, belong to orders with moderate to high diversification rates, which originated between 102 and 77 Mya, approximately between the Cenomanian and the mid-Campanian, in the Upper Cretaceous. According to constrained dating, the clades with the highest diversification rates appeared approximately between 95 and 75 Mya (mid-Cenomanian to mid Campanian, Upper Cretaceous) for stem orders, 80–30 Mya (mid-Santonian to latest Eocene) for crown orders, and 85–75 Mya (Santonian to latest Eocene) for nonnested crown clades. We emphasize that, although the differentiation and diversification of each order can be traced to an early time during angiosperm history, the timing and tempo of origination of the species contained in them cannot be provided by this study. The origin of the species contained in an order may have been concentrated at the onset of the diversification of the crown group, concentrated close to the present, or interspersed in-between in a variety of patterns. The data and analyses implemented in this study cannot discriminate among the numerous alternatives regarding within-lineage species-level diversification.

Absolute rates of diversification represent a starting point for a comprehensive understanding of the evolutionary mechanisms and causal factors that lead to the extraordinary extant angiosperm diversity. Investigation on the macroevolutionary dynamics of species accumulation in different angiosperm lineages are almost entirely lacking. Many questions still need to be investigated and resolved to provide integrative explanations to the abominable mystery. Did angiosperm species diversity originate at a gradual pace through time, or rapidly in a short burst sometime during the lineage's history? Are lineages with low or with high diversification rates more commonly characterized by time-homogeneous or time-heterogeneous processes? For some angiosperm lineages, including the earliest ones, the combination of an old age, a very small extant species diversity, and a well-documented fossil record provides circumstantial evidence in favor of a time-heterogeneous diversification process involving a rapid burst of species accumulation at the onset of the lineage's evolution. Explicit information about the timing and evolutionary mechanisms of angiosperm species origination would ideally require species-level dated phylogenies combined with an abundant fossil record and the ability to vary speciation and extinction rates across the tree and through time.

Previous studies (72; 53; 77) strongly imply that traits that conferred high diversification rates most likely evolved independently in different angiosperm lineages. We also believe it is unlikely that the outstanding angiosperm evolutionary success can be attributed to a single or a few particular intrinsic traits. Rather, it seems that many different traits, alone or in combination, placed at the right time under appropriate environmental conditions, allowed some groups within angiosperms to diversify spectacularly, and almost certainly at different times throughout their history.