Nested radiations and the pulse of angiosperm diversification

14 Pages·2015·1.11 MB·English

Research
Nested radiations and the pulse of angiosperm diversification:
increased diversification rates often follow whole genome
duplications
DavidC.Tank1,2,JonathanM.Eastman1,2,MatthewW.Pennell1,2,PamelaS.Soltis3,DouglasE.Soltis3,4,
CodyE.Hinchliff5,JosephW.Brown5,EmilyB.Sessa4andLukeJ.Harmon1,2
1DepartmentofBiologicalSciences,UniversityofIdaho,Moscow,ID83844,USA;2InstituteforBioinformaticsandEvolutionaryStudies,UniversityofIdaho,Moscow,ID83844,USA;
3FloridaMuseumofNaturalHistory,UniversityofFlorida,Gainesville,FL32611,USA;4DepartmentofBiology,UniversityofFlorida,Gainesville,FL32611,USA;5DepartmentofEcology
andEvolutionaryBiology,UniversityofMichigan,AnnArbor,MI48109,USA
Summary
Authorforcorrespondence: (cid:1) Ourgrowingunderstandingoftheplanttreeoflifeprovidesanovelopportunitytouncover
DavidC.Tank themajordriversofangiospermdiversity.
Tel:[phone] (cid:1) Usingatime-calibratedphylogeny,wecharacterizedhotandcoldspotsoflineagediversifi-
Email:[email protected]
cation across the angiosperm tree of life by modeling evolutionary diversification using
Received:10August2014
stepwise AIC (MEDUSA). We also tested the whole-genome duplication (WGD) radiation
Accepted:1May2015
lag-time model, which postulates that increases in diversification tend to lag behind estab-
lishedWGDevents.
NewPhytologist(2015)207:454–467 (cid:1) Diversification rateshavebeenincredibly heterogeneousthroughouttheevolutionaryhis-
doi:10.1111/nph.13491 toryofangiospermsandrevealapatternof‘nestedradiations’–increasesinnetdiversification
nestedwithinotherradiations.Thispatterninturngeneratesanegativerelationshipbetween
cladeageanddiversityacrossbothfamiliesandorders.Wesuggestthatstochasticallychang-
Keywords: angiospermdiversificationrates,
cladeageanddiversity,modeling ingdiversificationratesacrossthephylogenyexplainthesepatterns.Finally,wedemonstrate
evolutionarydiversificationusingstepwise significantstatisticalsupportfortheWGDradiationlag-timemodel.
AIC(MEDUSA),nestedradiations,poly- (cid:1) Across angiosperms, nested shifts in diversification led to an overall increasing rate of net
ploidy,whole-genomeduplicationradiation
diversificationanddecliningrelativeextinctionratesthroughtime.Thesediversificationshifts
lag-timemodel.
areonlyrarelyperfectlyassociatedwithWGDevents,butcommonlyfollowthemafteralag
period.
Thevaryingtempo of angiosperm diversificationhaslong fas-
Introduction
cinated biologists (Darwin, 1903). This mystery has only been
Withover 250000species (Crane etal.,1995;Juddetal.,2007; compounded as many of the recalcitrant branches of the angio-
and probably many more: e.g. Joppa etal., 2011), angiosperms spermtreeoflifehavebeenresolvedandacomprehensiveviewof
areamongthemostimpressiveradiationsofterrestrialorganisms. relationships among the major lineages of angiosperms has
In addition to their incredible species diversity, angiosperms are emerged(Soltisetal.,2011).Atthesametime,methodsforiden-
ecologically, functionally, and morphologically diverse, having tifying shifts in the rate of net diversification (speciation minus
risentodominanceinmostterrestrialenvironmentsina(geologi- extinction) have advanced considerably, allowing us to pinpoint
cally)shortperiodoftimesincetheirmajordiversificationinthe ‘hot’ and ‘cold’ spots on phylogenetic trees (e.g. where rates of
Cretaceous(Lidgard&Crane,1990;Craneetal.,1995;Dilcher, speciation and extinction have changed; Alfaro etal., 2009).
2001).Withinthismassiveradiationofspecies,itremainsamys- Despite decades of work characterizing angiosperm diversifica-
tery why some plant groups are much more diverse than others. tion from both paleontological and phylogenetic perspectives
Forexample,thecosmopolitanfamilyAsteraceae(whichincludes (Lidgard&Crane,1990;Sanderson &Donoghue, 1994;Crane
>23000species,includingdaisiesandsunflowers,inawidevari- etal.,1995;Dilcher,2001;Magall(cid:1)on&Sanderson,2001;Davies
ety of habitats) is <50 million yr old (Bremer & Gustafsson, & Barraclough, 2004; Magall(cid:1)on & Castillo, 2009; Smith etal.,
1997; Kim, 2005; Bell etal., 2010; Beaulieu etal., 2013), while 2011),westilldonothaveaclearideaofthedriversofdifferen-
its sister group, Calyceraceae, is restricted to South America and tial diversification across plant species. There have been many
comprisesonly60species(Stevens,2001-onwards;Bremeretal., hypothesesputforwardtoexplainthedramaticrisetodominance
2004). of angiosperms, and explanations are often based on ecological
454 NewPhytologist(2015)207:454–467 (cid:1)2015TheAuthors
www.newphytologist.com NewPhytologist(cid:1)2015NewPhytologistTrust
New
Phytologist Research 455
andphysiologicalinnovationsandthefunctionaltraitsassociated diversification rates across lineages generate this pattern (Stadler
withthese(Berendse&Scheffer,2009;Brodribb&Feild,2010; etal., 2014). Second, we show that increases in diversification
Labandeira, 2010; Feild etal., 2011; reviewed in Augusto etal., tend to lag behind major paleopolyploid events. These delayed
2014).Recently,withtheavailabilityofgenome-scaledatarepre- bursts may suggest that genome duplications promote, but are
sentinganincreasinglywidevarietyofangiospermlineages,anda not sufficient to cause, increased diversification. We suggest
stabilizing phylogenetic framework to investigate genome evolu- how such duplications might interact with other biotic and abi-
tion inacomparative framework,there hasbeenincreased inter- otic influences to lead to the nested pattern of radiations that
est in the role of paleopolyploidy – or ancient whole-genome we see across plantspecies.
duplication (WGD) – with respect to increased rates of lineage
diversification in angiosperms (Schranz etal., 2012; Vanneste
MaterialsandMethods
etal.,2014).
Recentgenomicinvestigationsindicatethatpolyploidyisubiq-
Methodologicaloverview
uitous among angiosperms, with evidence for an ancient WGD
event preceding the origin of the clade itself (Jiao etal., 2011; We present fully detailed explanations of our methods in the
Amborella Genome Project, 2013). Genomic data also suggest following sections. As a general overview, we first extended a
other major ancient WGD events in angiosperms (Van de Peer time-calibrated analysis of the Soltis etal. (2011) angiosperm
etal., 2009, 2010), including two WGDs that occurred in close phylogeny (Zanne etal., 2014) to incorporate phylogenetic
temporal succession early in eudicot evolution (Jiao etal., 2012) uncertainty. We then used the modeling evolutionary diversi-
withothereventsclosetotheoriginofmonocots,Poales,andSo- fication using stepwise AIC (MEDUSA) approach developed
lanales (Soltis etal., 2009). At least 50 independent ancient by Alfaro etal. (2009; extended by Pennell etal., 2014) to
WGDs are distributed across flowering plant phylogeny (Cui identify a set of clade-specific shifts in diversification rate
etal.,2006;Soltisetal.,2009;VandePeeretal.,2009,2010;M. across angiosperms. Finally, we placed nine well-characterized
S. Barker, pers. comm.). With this increased understanding of angiosperm WGDs on the phylogeny to test the WGD radia-
the ubiquitous nature of polyploidy in angiosperm histories and tion lag-time model (Schranz etal., 2012) of angiosperm
the phylogenetic location of these paleopolyploid events, it has diversification.
been recognized that a common pattern in the tree of life – spe-
cies-poor lineages subtending species-rich crown groups – tends
Time-calibratedphylogeny
to follow established WGD events in angiosperms (Soltis etal.,
2009, 2014b; Schranz etal., 2012). Schranz etal. (2012) erected We obtained a phylogenetic tree with branch lengths propor-
a formal hypothesis – the WGD radiation lag-time model – tional to time from a previous publication (Zanne etal., 2014).
based on observations of tree imbalance in several angiosperm Full details of tree construction are available in the original
lineages that have been associated with paleopolyploidy, which paper. Briefly, we compiled sequences for eight plastid genes
postulates that upticks in diversification rates tend to follow (atpB, matK, ndhF, psbBTNH, rbcL, rpoC2, rps16 and rps4)
WGD events, but only after lag times that may span millions of from Soltis etal. (2011), and used a by-gene partitioned maxi-
years (Schranz etal., 2012). However, there has been no formal mum likelihood analysis in RAXML v. 7.4.1 (Stamatakis, 2006;
quantitative test of this hypothesis across the angiosperm tree of Ott etal., 2007) to obtain a 639-taxon tree. Tree topology was
life. constrained according to Soltis etal. (2011) to ensure concor-
Our growing understanding of the structure of the plant tree dance with well-supported deep relationships among taxa
of life provides a rare opportunity to investigate angiosperm despite (relatively) limited genetic data. This ML tree was time-
diversificationtotrytouncoverthemajordriversthathaveledto scaled using 39 fossil calibrations and a penalized likelihood
thediversityofplantspeciesonEarth.Theconvergenceofawell- approach (see Zanne etal., 2014 for full details). The fossil cali-
sampled, well-resolved, and well-supported angiosperm phylog- brations used by Zanne etal. (2014) represent a reliable set of
enyandthemethodsfordetectingheterogeneityindiversification fossils that could be confidently identified and placed on the
rates from phylogenetic information – including topology, phylogeny, and, in addition to Zanne etal. (2014), these fossils
branchlengths,andcladerichnesses–inasinglecoherentframe- have all been used in other comprehensive large-scale dating
worksetsthestageforaddressingthequestion‘whatdrivesdiffer- analyses in plants (Bell etal., 2010; Smith etal., 2010; Beaulieu
entialdiversificationratesandtheresultingstrikingdisparitiesin etal., 2013). Rate smoothing was conducted by penalized likeli-
cladediversitiesacrossangiosperms?’ hood (TREEPL; Smith & O’Meara, 2012) using a smoothing
Here we use phylogenetic and taxonomic data to investigate parameter of 0.1 that was optimized on the maximum likeli-
angiosperm diversification dynamics. First, we characterize hot hood tree. For each calibration, both minimum and maximum
and cold spots of diversification across the angiosperm tree of age constraints were applied, where minimum age constraints
life. We identify a heterogeneous pattern of rapid radiations corresponded to the age of the fossil used in previous analyses
nested within other rapid radiations, which we refer to as (Bell etal., 2010; Smith etal., 2010; Beaulieu etal., 2013) and
‘nested radiations’ – similar to the ‘repeated radiations’ observa- maximum age constraints were calculated from the upper
tion of Soltis etal. (2004), but here with formal diversification 97.5% of the lognormal distribution with means and standard
rate analyses – and suggest that stochastically changing deviations following the lognormal priors used for the Bayesian
(cid:1)2015TheAuthors NewPhytologist(2015)207:454–467
NewPhytologist(cid:1)2015NewPhytologistTrust www.newphytologist.com
New
456 Research Phytologist
divergence time estimates of Bell etal. (2010), Smith etal. computing the average and standard deviation of the mean rates
(2010), and Beaulieu etal. (2013; see Zanne etal., 2014 supple- from each tree, as well as the overall range of the parameter
mentary Table 2 for full details). In addition to the fossil cali- estimates across all the trees. When resolution is lost at the fam-
brations, Zanne etal. (2014) constrained the root node with a ily level (as a result of the incorporation of unsampled diver-
minimum age of 301 million yr (Myr) and a maximum of sity), the average parameter values represent only those lineages
366Myr following the results of Smith etal. (2010) and recom- with resolution, resulting in an increasingly smaller number of
mendations of Clarke etal. (2011). To account for uncertainty lineages contributing to the calculation of mean rates towards
in divergence time estimation, we repeated this analysis pipeline the present.
for a set of 1024 bootstrapped data sets. For each of these repli-
cates, we resampled the genetic sequence data with replacement
Cladeage–diversityrelationships
to obtain new alignments with the same length as the original,
and then repeated all of the analyses described above to obtain a We tested for a relationship between clade age and diversity at
chronogram with branch units in millions of years. We applied the level of plant families and orders. To do this, we correlated
all the diversification rate analyses described in the following the species diversity of each family (n=325) and each order
section across thisfull set of 1024bootstrappedchronograms. (n=66)withitsstemageusinglinearregression.
Identifyingdiversificationshifts Whole-genomeduplicationsandcorrelationwith
diversificationrates
We used MEDUSA (Alfaro etal., 2009; Pennell etal., 2014) to
identify shifts in diversification rates along branches in the seed Finally, we explored the relationship between WGDs and shifts
plant phylogeny. After initially fitting a constant rate birth– indiversificationrates.Todothis,weplacedninewell-character-
death model of diversification to the phylogenetic tree, izedWGDsatthefamilylevelandaboveonthetreeandaskedif
MEDUSA uses a step-wise addition algorithm to infer phyloge- thereisacorrespondencebetweenpolyploidizationandincreases
netically local shifts in the rates of two diversification parame- in net diversification rate. We selected these nine WGDs based
ters: net diversification (r=k(cid:3)l) and relative extinction onboth thestrength ofevidencein support of aWGDatapar-
(e=l/k), where k is the speciation rate (birth) and l is the rate ticular node and our ability to precisely place a WGD on our
of extinction (death). Rate shifts are retained if including the phylogenygiventhesampling.Becausewecollapsedcladesatthe
shift substantially improves the sample-size corrected Akaike familyleveltoincorporateunsampleddiversityinourdiversifica-
Information Criterion (AIC) score (AICc; Burnham & Ander- tion rate analyses, we cannot place any WGD events within a
son, 2002). The version of MEDUSA we used for this analysis family. In addition, because the hypothesized phylogenetic loca-
has been improved from the original implementation in three tionsofmanyWGDsinferredfromgenomicdataareoftenbased
ways: it (1) considers a mixture of pure-birth and birth–death on very limited taxonomic sampling, it is difficult to precisely
processes in shifted lineages; (2) uses an AIC threshold that is placetheseeventsonaphylogeny.
corrected for tree size; and (3) considers both forward and Table1showsthenineWGDeventsthatwereusedtotestfor
backward selection in model choice (see Pennell etal., 2014 for a correlation with increased rates of diversification. In all cases,
more details). theseeventswereidentifiedthroughcomparativegenomicanaly-
For this MEDUSA analysis, we collapsed our 639-taxon phy- ses. The first two represent established WGDs that occurred in
logeny into anexemplartree,with onerepresentative foreachof theancestorofallangiosperms(theeevent;Jiaoetal.,2011)and
the325sampledseedplantfamilies.BecauseMEDUSArequires the subsequent radiation of the core eudicots (the c event; Jiao
that all missing species be assigned to a tip clade in the full tree, etal., 2012). The third established paleopolyploidization repre-
we mapped species richnesses to each familial exemplar using a sentsanancientgenomeduplicationinmonocots.Basedoninte-
comprehensive systematic resource, The Plant List (2010). We grated syntenic and phylogenomic analyses, Jiao etal. (2014)
wereabletoaccountforalargeproportion(0.949)ofknownspe- revealed a WGD shared by all monocots sampled to date (the s
ciesrichnessinseedplants.Inferencesofrateheterogeneityinthe event).Itisimportanttonotethatthiseventisseparatefromthe
diversificationprocessofseedplantsaredrawnfromthecompila- previously characterized q and r events in monocots (Paterson
tion ofMEDUSAanalyses acrossour fullbootstrappeddistribu- etal., 2012). However, because no noncommelinid monocot
tion of time-calibrated trees (n=1024; Supporting Information genomeshavebeensampledthusfar,thephylogeneticplacement
TableS1). of the s event is more uncertain. Jiao etal. (2014) hypothesized
To summarize the time-course of speciation and extinction thatsmayrepresentaWGDinoneoftheearlydivergingmono-
from our MEDUSA analysis, we divided each tree into 1-Myr cot lineages, but given that sampling to date has been limited
time intervals. For each time interval, we calculated the mini- only to the Commelinidae (i.e. palms, bananas, and grasses), we
mum, maximum, and mean of all branch-specific rate estimates considered several alternative placements for this event that
of net diversification (r) and relative extinction (e) for all ranged from the commelinid clade, which is the most recent
branches occurring in that interval. This results in a time-series place this could have occurred given the available data, to the
of rate estimates for each tree. We then repeated this procedure entire monocot clade, including intermediate alternative place-
over all of the bootstrap trees, and summarized the results by ments along the backbone of the clade (Table1). In their
NewPhytologist(2015)207:454–467 (cid:1)2015TheAuthors
www.newphytologist.com NewPhytologist(cid:1)2015NewPhytologistTrust
New
Phytologist Research 457
Table1 Whole-genomeduplications(WGDs)investigatedforcorrespondencewithincreaseddiversificationrates
WGDdescription Reference Nearestrincrease n-dist Myr
1 Angiospermae(e) Jiaoetal.(2011) Mesangiospermae 3 49.2
2 Gunneridae(c) Jiaoetal.(2012) MRCAof(Superrosidae+Superasteridae) 1 0.6
3 Commelinidae(s) Jiaoetal.(2014) MRCAof(Arecaceae)and 1 2.3
(Commelinales+Zingiberales)
3a Monocotyledonae(s) Jiaoetal.(2014) MRCAof(Commelinidae+Asparagales) 4 39.4
3b Nartheciidae(s) Jiaoetal.(2014) MRCAof(Commelinidae+Asparagales) 3 26.5
3c Petrosaviidae+Dioscoreales(s) Jiaoetal.(2014) MRCAof(Commelinidae+Asparagales) 2 5.7
3d Petrosaviidae(s) Jiaoetal.(2014) MRCAof(Commelinidae+Asparagales) 1 3.0
3e Commelinidae+Asparagales(s) Jiaoetal.(2014) MRCAof(Commelinidae+Asparagales) 0 0.0
4 Poaceae(q,r) Patersonetal.(2012); NA NA NA
Jiaoetal.(2014)
5 Magnoliales+Laurales Cuietal.(2006) Lauraceae 4 26.8
6 Ranunculales Cuietal.(2006); None NA NA
Pabo(cid:1)n-Moraetal.(2013)
6a Papaveraceae+Ranunclaceae Pabo(cid:1)n-Moraetal.(2013) None NA NA
7 Asteraceae Barkeretal.(2008) Asteraceae 0 0.0
8 Brassicaceae(a) Haudryetal.(2013); NA NA NA
Kagaleetal.(2014)
9 Brassicaceae+Limnanthaceae(b) Barkeretal.(2009) MRCAof(Capparaceae+Brassicaceae) 3 35.9
9a CoreBrassicales Barkeretal.(2009) MRCAof(Capparaceae+Brassicaceae) 1 19.9
[Resedaceae+Brassicaceae](b)
9b Brassicaceae+Bataceae(b) Schranzetal.(2011) MRCAof(Capparaceae+Brassicaceae) 2 28.8
WGDdescription,locationofWGD,withnamesfollowingSoltisetal.(2011);reference,citationforeachWGD;nearestrincrease,locationofthenearest
increaseinnetdiversificationthatwasidentifiedin>75%oftrees(r);n-dist,thenumberofnodesbetweenWGDsandthenearesttip-wardincreaseinnet
diversification;Myr,thelagtime(inmillionsofyears)betweenWGDandshiftedcladecalculatedusingdivergencetimeestimatesfromthemaximumlikeli-
hoodtree.Rows3(a–e),6(a),and9(a,b),representalternativeplacementsoftheseWGDs(seetext).MRCA,mostrecentcommonancestor;NA,notappli-
cable.
analyses,Jiaoetal.(2014)alsodisambiguatetheplacementofthe duplicationoccurredbeforeorafterthedivergenceofEupteleaceae
q and r events, confidently placing both of these WGDs along (Pab(cid:1)on-Moraetal.,2013).Therefore,fortheranunculiidWGD,
thelineageleadingtograsses(Table1). we investigated an alternative placement excluding Eupteleaceae,
Inapioneeringcomparativegenomicstudyofpaleopolyploidy, thesistergrouptotherestoftheclade.
Cui etal. (2006) identified several WGD events in angiosperms, TheseventhWGDthatweconsideredisthewell-studiedpoly-
includingasharedWGDbetweenMagnolialesandLaurales,and ploidization event that characterizes the sunflowers (Asteraceae;
an additional duplication in Ranunculales. While these WGDs Barker etal., 2008). Finally, we also considered the a and b
have been widely cited (Soltis etal., 2009; Van de Peer etal., WGDs that have been extensively studied in Brassicales (Bowers
2009),giventhelimitedsamplinginthisstudy(i.e.asingletaxon etal., 2003; Barker etal., 2009; Schranz etal., 2011; Haudry
representing each order), generalizing these duplications to the etal.,2013;Kagaleetal.,2014).Whilethephylogeneticlocation
ordinal level is tentative. Using preliminary data from the 1000 oftheaeventhasbeenconfidentlyplacedonthelineageleading
plants(1KP)transcriptomeproject,wewereabletoconfirmboth toBrassicaceae(Haudryetal.,2013;Kagaleetal.,2014),because
of these events (P. S. Soltis & D. E. Soltis, unpublished data). oflimitedgenomicsamplinginBrassicales, theplacementofthe
FromtranscriptomesforfourspeciesofMagnolialesandsixspe- b WGD is more problematic. Based on the sampling from
cies of Laurales, initial plots of the divergence of duplicate gene Barker etal. (2009; Brassicaceae, Cleomaceae, and Caricaceae),
pairsintermsofsubstitutionspersynonymoussiteperyear(K), thedeepestplacementofthiseventisfollowingthedivergenceof
s
indicate a shared WGD that is unique to these orders. Likewise, Caricaceae (Brassicaceae+Limnanthaceae in Table1). However,
K analyses from seven species of Ranunculales revealed a shared Barker etal. (2009) hypothesize that this duplication probably
s
WGDassuggestedbyCuietal.(2006).However,becausediploi- representsapaleopolyploidizationofthecoreBrassicales(Reseda-
dization can obscure these patterns – especially at deep phyloge- ceae+Brassicaceae in Table1), and Schranz etal. (2011) suggest
netic levels – these data are part of an ongoing, more detailed thatthebWGDoccurredfollowingthedivergenceofLimnanth-
studythatwillintegratethesecomparativegenomicanalyseswith aceae(Brassicaceae+BataceaeinTable1).Therefore,wechoseto
syntenicandphylogenomictechniques(e.g.asinJiaoetal.,2014). include these three alternatives as possible placements for the b
Nevertheless, here we include a WGD eventsharedby Magnoli- WGD(Table1).
alesandLaurales,aswellasaRanunculalesWGD(Table1).For GiventheseWGDs,wefirstaskediftherewasaperfectcorre-
Ranunculales, additional evidence from phylogenetic analyses of spondence between inferred increases in net diversification and
MADS-boxgenesinRanunculalessupportsaduplicationearlyin polyploidization events. To accomplish this, we first determined
thediversificationofRanunculales,butitisnotclearwhetherthis how many of our nine identified WGD events corresponded
(cid:1)2015TheAuthors NewPhytologist(2015)207:454–467
NewPhytologist(cid:1)2015NewPhytologistTrust www.newphytologist.com
New
458 Research Phytologist
exactlywithoneofthe15inferredincreasesinnetdiversification estimatestoaccountforphylogenetic(temporal)uncertaintywith
rate from our MEDUSA analyses (Table2). We considered sev- respect to penalized likelihood analyses. Stem ages for angio-
eral alternative placements for three of the WGDs events spermfamilies(andassociated95%confidenceintervals)retained
(Table1); we tested all 36 (69293) combinations of these forourdiversificationanalysesareshowninTableS2.
alternative placements. For each set of WGD placements, we
counted the number of corresponding diversification rate shifts,
Variationindiversificationratesamongclades
and we compared this number with how many correspondences
one wouldexpectunder asimple null modelwhere nine WGDs Our analyses identified a total of 41 shifts in the maximum
wereplacedrandomlyonbranchesinthetree. likelihood estimate (MLE) tree and 142 total unique shifts in
WethentestedtheWGDradiationlag-timehypothesis,which diversification rate across the set of bootstrapped trees (Table
predictsadelayedelevatednetdiversificationratefollowingpoly- S1). However, incorporating uncertainty in divergence times
ploidization(Schranz etal.,2012). Inthiscase,ourhypothesisis reduces this number to 27 major shifts that occur in at least
that polyploidization events will elevate net diversification rates 75% of the bootstrap set (Table2; Fig.1). Shifts were distrib-
after some delay (measured here as nodal distance). We chose uted throughout the tree, revealing heterogeneous patterns of
nodal distance (as opposed to absolute time) to evaluate this diversification across angiosperms, but with most shift points
hypothesis, because there is no expectation of exactly how long clustered in the denser parts of the tree between 50 and 125
thelagtimeshouldbe.TheverbalmodelofSchranzetal.(2012) million yrago(Ma)(Fig.2a,b).Theseshifts often showanested
recognizes that the delayed shift in diversification rate may take pattern of rate shifts within a clade followed by further shifts
many millions of years. In addition, the precise temporal loca- within asubcladeof thatoriginalclade.
tions of both the identified diversification rate shifts and the Our calculated diversification rates through time show high
WGDsareunclear.Wevisualizebothoftheseasbeingatanode, variability among branches and across the bootstrap distribution
wheninrealitytheyhappenedsomewherealongthestemlineage. oftrees.Ingeneral,werevealedatrendofincreasingnetdiversifi-
NeitherMEDUSAnortheapproachesusedtoidentifythephylo- cation rates (r) and high, but decreasing, turnover (high e) rates
geneticpositionofWGDsprovidemeaningfultimeestimatesfor throughtime(Fig.2;Tables2,S1).Becausecladeswerecollapsed
where these events happened. Therefore, we chose to use nodal at the family level to incorporate unsampled diversity, there is a
distancewithvariouscut-offsthatrepresentareasonableexpecta- markedlevelingoffofthemeanrateestimatesforbothnetdiver-
tionforthetime-laghypothesis.ThreeoftheWGDevents(Poa- sification and relative extinction using this summary approach
ceae, Brassicaceae, and Asteraceae; Table1) are tip lineages. In (Fig.2c,d).
thesecases,wecannottesttheWGDradiationlag-timehypothe-
sis,becausewehavenoresolutioninthetreewithintheselineages
Cladeage–diversityrelationships
(i.e. the tree was collapsed to the family level for our diversifica-
tion rate analyses); we excluded these tip lineages from the lag- Consistent with previous studies of angiosperm diversification
timetest.ForeachoftheothersixWGDevents,wedetermined (Magall(cid:1)on & Sanderson, 2001; Magall(cid:1)on & Castillo, 2009), we
whetherornotanyofthe15identifiedshiftsinelevateddiversifi- foundasignificantnegativerelationshipbetweenageandthenat-
cation rates (Table2) followed within a set number of nodes on ural logarithm of species richness considered at the familial level
the tree (we used three nodes as our cut-off, but results were (R2=0.016; P=0.0227; Fig. S1A) and at the ordinal level
robusttocut-offsofonethoughfour,whichisthemaximumdis- (R2=0.144;P=0.0017;Fig.S1B).Stem-groupagesforthe325
tanceofanyrateshiftfromaWGDevent).Weusedthenumber angiosperm families summarized across our distribution of trees
ofWGDsthatwerefollowedbyarateshiftasateststatistic,and areprovidedinTableS2.
compared this to a null distribution generated by again placing
six WGDs randomly on branches in the tree. As before, we
WGDsanddiversificationrate
repeatedthistestforall36setsofpotentialWGDplacements.
We first tested for a perfect correspondence between diversifica-
tion-rate shifted lineages and polyploidization events. We found
Results
thateitheroneortwoofnineinvestigatedpolyploidizationnodes
were perfectly associated with diversification upticks (Fig.3;
Dataavailabilityandanalysispipelines
Table1; Asteraceae and Commelinidae+Asparagales, when this
Alldata setsandR-scripts necessary toperformtheanalyses pre- putativeplacementwasselectedtorepresentofthesevent(3ein
sentedhereareavailableathttps://github.com/harmonlab/angio- Table1)). We can only reject the null hypotheses in the six sce-
pulse. narios(outof36)withtwoexactmatches.
We then tested the WGD radiation lag-time hypothesis by
testing for a delayed correspondence between diversification-rate
Time-calibratedphylogeny
shiftedlineagesandpolyploidization(excludingthreetiplineages
Our extension of the divergence time estimates of Zanne etal. where we would be unable to detect a delayed upturn). Using a
(2014) resulted in a bootstrap set of 1024 time-calibrated trees. cut-off ofthreenodes,wefound thateither threeorfour(outof
This bootstrap set was used as a distribution of divergence time a possible six) polyploidization nodes show delayed
NewPhytologist(2015)207:454–467 (cid:1)2015TheAuthors
www.newphytologist.com NewPhytologist(cid:1)2015NewPhytologistTrust
New
Phytologist Research 459
Table2 MEDUSA(modelingevolutionarydiversificationusingstepwiseAIC)estimatesforasetofprimaryshiftsindiversificationfindingsupportinatleast
75%ofbootstrapreplicates
Description Label Support Ma r e Dr De
MRCAof(Capparaceae+Brassicaceae) 26 0.94 34.7 0.2187 0.00 0.1614 (cid:3)0.31
Cactaceae 27 1.00 29.1 0.2626 NA 0.1359 NA
Asteraceae 25 1.00 48.8 0.2101 NA 0.1292 NA
MRCAof(Caryophyllaceae)and(Stegnospermataceae) 23 1.00 62.2 0.1267 0.00 0.0910 (cid:3)0.03
MRCAof(Gentianales+Solanales)and(Boraginaceae+Lamiales) 19 0.90 76.3 0.1287 0.52 0.0821 0.61
Fabaceae 21 1.00 69.3 0.1452 NA 0.0809 NA
Piperaceae 20 0.98 74.8 0.1012 NA 0.0534 NA
MRCAof(Araliaceae)and(Myodocarpaceae+Apiaceae) 18 0.99 76.5 0.0996 0.00 0.0469 0.00
MRCAof(Vochysiaceae+Myrtaceae)and(Melastomataceae) 15 0.76 88.2 0.1102 0.00 0.0451 (cid:3)0.44
Annonaceae 17 0.99 81.2 0.0927 NA 0.0451 NA
Lauraceae 13 0.98 93.0 0.0857 NA 0.0379 NA
Mesangiospermae 1 1.00 218.8 0.0473 0.22 0.0349 (cid:3)0.66
MRCAof(Superrosidae+Superasteridae) 7 0.98 119.2 0.0645 0.54 0.0174 0.35
MRCAof(Commelinidae+Asparagales) 4 0.92 135.6 0.0626 0.94 0.0154 0.71
MRCAof(Arecaceae)and(Commelinales+Zingiberales) 5 0.92 123.8 0.0629 0.00 0.0003 (cid:3)0.94
MRCAof(Campanulidae+Lamiidae) 8 1.00 106.2 0.0612 0.01 (cid:3)0.0030 (cid:3)0.54
MRCAof(Paeoniaceae)and(Altingiaceae) 9 0.98 103.6 0.0377 0.00 (cid:3)0.0277 (cid:3)0.50
MRCAof(Rhabdodendraceae)and(Simmondsiaceae) 16 0.83 87.5 0.0293 0.00 (cid:3)0.0342 (cid:3)0.65
Ceratophyllaceae 2 0.95 186.1 0.0075 NA (cid:3)0.0401 NA
Huerteales 14 0.78 90.5 0.0223 0.00 (cid:3)0.0409 (cid:3)0.67
Acoraceae 3 0.97 172.4 0.0040 NA (cid:3)0.0435 NA
Blandfordiaceae 6 1.00 119.6 0.0116 NA (cid:3)0.0505 NA
MRCAof(Aphloiaceae+Strasburgeriaceae) 12 1.00 93.9 0.0074 0.00 (cid:3)0.0538 (cid:3)0.63
Berberidopsidales 10 1.00 97.9 0.0043 0.00 (cid:3)0.0599 (cid:3)0.55
Curtisiaceae 11 0.96 96.1 0.0000 NA (cid:3)0.0645 NA
MRCAof(Montiniaceae)and(Sphenocleaceae+Hydroleaceae) 22 1.00 63.1 0.0324 0.00 (cid:3)0.0984 (cid:3)0.68
Plocospermataceae 24 0.92 61.5 0.0000 NA (cid:3)0.1288 NA
Whereindicated,valuesareaveragedacrossbootstrap-replicatetrees.
Description,topologicalpositionoftheshiftedlineage,withcladenames(initalics)followingCantinoetal.(2007)andSoltisetal.(2011),andfamilyand
ordernamesfollowingAPGIII(TheAngiospermPhylogenyGroup,2009);label,temporalorderingofeachshift,beginningwiththeoldest(i.e.
Mesangiospermae);support,fractionoftreesinwhichshiftreceivedstatisticalsupport;Ma,meanestimateofageofshift(millionsofyearsago);r,maxi-
mumlikelihoodestimateofnetdiversification;e,maximumlikelihoodestimateofrelativeextinction;D,meanmagnitudeofchangeinnetdiversification
r
relativetoimmediateancestor(seeFig.1aforfurtherdetails);De,meanmagnitudeofchangeinrelativeextinctionrelativetoimmediateancestorrelative
toimmediateancestor(seeFig.1aforfurtherdetails).Becausediversificationrateanalyseswerecollapsedtothefamilylevel,relativeextinction(eandDe)
cannotbecalculatedandisnotedasNA(notapplicable)forthoseclades.MRCA,mostrecentcommonancestor.
diversification upticks (Fig.3; Table1; Angiospermae, Later work used a variety of methods to locate and characterize
Gunneridae, the monocot s event in all positions except 3b, and these shifts (Magall(cid:1)on & Sanderson, 2001; Davies & Barrac-
all placements of the Brassicales b event). In all 36 cases we can lough,2004;Magall(cid:1)on&Castillo,2009;Smithetal.,2011;Fiz-
rejectthenullhypothesis(P<0.001–0.011),andusingthealter- Palaciosetal.,2011).Incontrasttothesestudies,theanalyseswe
nate cut-offs of one, two, orfournodes,we could reject thenull performed here utilized all aspects of the phylogenetic data –
hypothesis of a random association in 20, 30, or all 36 cases, includingthetopology,branchlengths,andcladerichnesses–in
respectively, and thus support a delayed association between asinglecoherentframeworkbasedonbirth–deathmodels.
polyploidy and an uptick in diversification rate that ranges in Heterogeneity in diversification rates may occur either tem-
time from 0.6Ma (Gunneridae WGD and porally (e.g. speciation rates slowing through time) or among
Superasteridae+Superrosidae shift) to 49.2 Ma (Angiospermae lineages. Our analyses focus on the latter, as MEDUSA is a
WGDandMesangiospermaeshift;Fig.3;Table1). clade-based approach used to test for shifts in diversification
rates between lineages (Alfaro etal., 2009). An assumption of
temporal-based methods is that changes in diversification rates
Discussion
occur homogeneously across the tree. However, for a large
group such as angiosperms, this assumption is almost certainly
Heterogeneousdiversification
violated: previous studies (Magall(cid:1)on & Sanderson, 2001;
We found a striking amount of heterogeneity in diversification Davies & Barraclough, 2004; Magall(cid:1)on & Castillo, 2009;
rates across angiosperms (Table2; Figs1, 2). Previous work has Smith etal., 2011) have found substantial differences in diver-
suggestedthattheoverallincreaseddiversificationratesinangio- sification rates among angiosperm clades. Furthermore, our
spermsareassociatedwithrateshiftsthatoccurredaftertheinitial level of sampling at the species level is both low (639 exemplar
diversification of angiosperms (Sanderson & Donoghue, 1994). taxa out of >250000 known angiosperms) and highly
(cid:1)2015TheAuthors NewPhytologist(2015)207:454–467
NewPhytologist(cid:1)2015NewPhytologistTrust www.newphytologist.com
New
460 Research Phytologist
(a)(b) eSummaryofMEDUSA(modelingevolutionarydiversificationusingstepwiseAIC)parameterestimatesandshiftsfornetdiversification(r;a)andrelativeextinction(;b),plottedontherate-Fig.1smoothedmaximumlikelihoodestimate(MLE)ofthephylogeny(seetext).Collapsedcladesaredrawnastrianglesthatareproportionaltofamily-levelspeciesrichnesses.Colorsofbrancheseecorrespondtoparameterestimates(seekeysforrand)fortheMLE;colorswithincirclesreflectmagnitudesofchangeinr(a)or(b)fromtheimmediatephylogeneticbackgroundtotheshiftedeDDlineage,averagedacrossthedistributionofbootstrapreplicates.Boldercolorsareassociatedwithmoreextremevaluesofrandorwithgreatermagnitudesofshifts(and).In(b),graywithinercirclesisusedwherea(unresolved)familyisestimatedtohaveindependentratesofdiversification.Estimatingrelativeextinctionisproblematicinthesecircumstancesbecauseofthelackofresolutioninthesubtree(Raboskyetal.,2007).ThesizeofcirclesindicatesshiftsupportastherelativefrequencyoftreesinthedistributionforwhichMEDUSArecoversagivenshiftedlineage(see‘support’key).Forclarity,supportforshiftsatunresolvedcladesisindicatedattheoutermostextentsofthoseclades(e.g.Cactaceae).Majorlineagesarelabeledinternally,andfamilieswithsupportforashiftingreaterthan5%ofanalysesacrossthedistributionoftreesareenlargedastiplabels.Thesesamefamilieshaverepresentativeicons,foundeitherin(a)orin(b).Allothertiplabelsareminimized.CircleswithnumericlabelscorrespondtothosefoundinFig.2andTable2.
NewPhytologist(2015)207:454–467 (cid:1)2015TheAuthors
www.newphytologist.com NewPhytologist(cid:1)2015NewPhytologistTrust
New
Phytologist Research 461
(a)
(b)
(c)
(d)
nonrandom so that inferences based on temporal patterns etal., 2009) by plotting the timing of rate shifts calculated
using available methodologies are likely to be misleading from 1-Myr time intervals across the trees (Fig.2). The major-
(Cusimano & Renner, 2010). ity of shifts in diversification rate occur in the interval between
Still, we attempted to summarize the temporal pattern of c. 125 and 50 Ma (Fig.2a,b), even when we accounted for
rate shifts inferred using the MEDUSA algorithm (Alfaro some of the uncertainty in estimating divergence times
(cid:1)2015TheAuthors NewPhytologist(2015)207:454–467
NewPhytologist(cid:1)2015NewPhytologistTrust www.newphytologist.com
New
462 Research Phytologist
Fig.2 SummaryofresultsfromMEDUSA(modelingevolutionarydiversificationusingstepwiseAIC),focusingon27shiftsthatarewellsupportedacrossa
distributionofbootstrapreplicates(i.e.foundin>75%ofbootstrapreplicates).(a)Meantimingandnatureofprimaryshifts:redpointsdepictmagnitude
ofchangeinrelativeextinctioncomparingashiftedlineage(es)toitsimmediateancestor(ea)suchthatDe=es(cid:3)ea.Similarly,blacksegmentsindicate
ancestor–descendantchangesinnetdiversification(D =r (cid:3)r).NumericlabelscorrespondtothoseinTable2andareorderedthroughtime;forexample,
r s a
thefirstshift(atc.219millionyrago(Ma),forMesangiospermae;Table2)showsamoderateincreaseinnetdiversificationandaprecipitousdropin
relativeextinctionincomparisontotheimmediatebackgroundrates.(b)Uncertaintyintimingofprimaryshifts,compiledfromMEDUSAanalysesforeach
treeinthedistributionofbootstrapreplicates.Temporalbinsof1millionyr(Myr)wereusedtoassesstheshiftdensitythroughtime(seetext).Broader
peaksarethoseforwhichthereismoreuncertaintyintiming(e.g.themesangiospermshiftbetweenc.230and200Ma).Densitiesmayexceed1if
multipleshiftsoverlapintiming.(c,d)Trendsthroughtimeinestimatedmeans(black)andranges(gray,wheredarkgrayrepresentsuncertaintyin
estimatedmeansacrossthebootstrappeddistributionoftrees,andlightgrayrepresentstheabsoluterangesacrossalltrees)areshownfornet
diversification(c)andrelativeextinction(d).Parameterestimatesforallextantlineagesinagiven1-Myrbinoftime,andacrossallbootstrapreplicates,
wereusedtogeneratepanels(c)and(d).Asaresultoftheinabilitytoestimaterelativeextinctionincladesthatlackresolution(Raboskyetal.2007),
estimatesassociatedwithshiftsinunresolvedcladeswereexcludedfrom(d)and(a).
(Fig.2b). This is perhaps unsurprising as this interval coincides period coincides with a nonrandom association of WGD
with the rise of many of the major lineages of angiosperms, events (Vanneste etal., 2014).
especially among the eudicots. This period also encompasses Weobservetwogeneraltemporaltrendsinangiospermdiversi-
the Cretaceous–Paleogene (K-Pg) extinction event, which pa- fication.TheoriginofMesangiospermaemarksthebeginningofa
leobotanical studies suggest led to the extinction of up to c. trend toward gradually increasing rates of net diversification
60% of plant species in some regions (Wilf & Johnson, coupledwithatrendofdecreasingrelativeextinction(Fig.2c,d).
2004), and which provoked major changes in regional floras Following this, we further observe a trend toward increasing
(McElwain & Punyasena, 2007). Furthermore, this time among-lineage heterogeneity in both net diversification and
Fig.3 Phylogeneticlocationofwhole-genomeduplication(WGD)eventsandthecorrespondingincreasesinnetdiversificationthatsupporttheWGD
radiationlag-timemodel(Schranzetal.,2012).Collapsedcladesaredrawnastrianglesproportionaltofamily-levelspeciesrichnesses.Greenstarsindicate
WGDevents,andcoloredcirclesshowthelocationofincreasesinnetdiversification.Colorswithincirclesreflectmagnitudesofchangeinr(net
diversification)fromtheimmediatephylogeneticbackgroundtotheshiftedlineage,averagedacrossthedistributionofbootstrapreplicates;boldercolors
areassociatedwithgreatermagnitudesofshifts(D).Thesizeofcirclesindicatesshiftsupportastherelativefrequencyoftreesinthedistributionforwhich
r
MEDUSA(modelingevolutionarydiversificationusingstepwiseAIC)recoversagivenshiftedlineage(see‘support’key).Forclarity,supportforshiftsat
unresolvedcladesisindicatedattheoutermostextentsofthoseclades(e.g.Asteraceae).BothWGDsandassociateddiversificationrateshiftsare
numberedaccordingtoTable1,andmajorlineagesarelabeled.
NewPhytologist(2015)207:454–467 (cid:1)2015TheAuthors
www.newphytologist.com NewPhytologist(cid:1)2015NewPhytologistTrust
New
Phytologist Research 463
relative extinction. This large variance suggests that methods which is in turn determined by an odd mix of biology and psy-
implicitly assuming rate homogeneity across the tree may be chology. Stadler etal.’s (2014) simulations clearly demonstrate
inappropriate for this data set (but see Morlon etal., 2011). that it is possible to generate inverse relationships between clade
Many hypotheses have been put forward to explain the increase age and richnessunder models without ecological limits. Impor-
in diversification rates among angiosperm lineages (Stebbins, tantly, weak or negative age–diversity relationships are expected
1970, 1971, 1974, 1981; reviewed in Gorelick, 2001; Davies & under‘nestedradiation’conditionsasrevealedbyourdiversifica-
Barraclough,2004;Vamosi&Vamosi,2011;reviewedinAugu- tionrateanalyses.
sto etal., 2014), and we do not attempt to distinguish between Our arguments against an ‘ecological limits’ interpretation of
thesehere,buthopethatotherresearcherswillbeabletouseour theage–diversityrelationshipdonotmeanecologicalinteractions
resultstotestthesehypothesesandgeneratenewones. have been inconsequential to the diversification of angiosperms.
Modeling approaches that consider diversification rates as a Indeed, there is some evidence from the fossil record (Alroy,
functionoftime(Pybus&Harvey,2000;Rabosky,2006;Rabo- 2008,2010),andincreasinglyfromreconstructedmolecularphy-
sky & Lovette, 2008; Morlon etal., 2010, 2011; Stadler, 2011) logenies (Phillimore & Price, 2008; Rabosky & Lovette, 2008;
include some that can model variation both across clades and Rabosky&Glor,2010),thatdiversity-dependentspeciationand
through time (Morlon etal., 2011; Rabosky, 2014). Such extinctionhaveplayedaprominentroleinthegenerationofbio-
approaches (reviewed in Pennell & Harmon, 2013; Pyron & diversity. However, it is important to note that the presence or
Burbrink, 2013; Morlon, 2014) represent an interesting next- absence of ecological limits cannot be reliably judged from the
step analysis for angiosperm radiations, especially as more fully signoftherelationshipbetweenlineageageandrichness(Stadler
resolvedtime-treesemerge. etal.,2014).
Negativeage–diversityrelationshipacrossplantclades WGDsasadriverofnestedradiations
Across angiosperm clades, we found a significant negative rela- WGDshavelongbeenrecognizedasimportantdriversofspecia-
tionship between clade age and the natural logarithm of species tion in plants (Clausen etal., 1945; Stebbins, 1947, 1950), but
richness(Fig.S1;seealsoTableS2).Thisresultconfirmsthefind- studies explicitly linking rates of diversification and polyploidy
ingsofpreviousstudies(Magall(cid:1)on&Sanderson,2001;Magall(cid:1)on havebeenrelativelyfew(Mayroseetal.,2011;Soltisetal.,2014a;
& Castillo, 2009), but is notable given our more complete sam- see Zhan etal., 2014 for an example in fish), and these have
pling of angiosperm diversity included here and the updated focused primarily on the evolutionary consequences of recent
divergence time analyses. A negative or nonexistent correlation polyploid formation, with contrasting results (see Soltis etal.,
betweencladeageandcladerichnesshasalsobeenreportedfora 2014a for a review of these analyses). With the rise of genome-
number of other groups (e.g. avian tribes: Ricklefs, 2006; major scale sequencing studies and statistical methods for identifying
squamateclades:Ricklefsetal.,2007;snakes:Pyron&Burbrink, WGD events across angiosperms, we are gaining confidence in
2012),andonerecentanalysisprovidesevidencethatthismaybe theplacementofanumberofancestralgenomeduplications(e.g.
acommonfeatureacrossthetreeoflife(Raboskyetal.,2012). seed plants and angiosperms: Jiao etal., 2011; monocots: Jiao
A number of recent papers (Rabosky, 2009a,b, 2010, 2012; etal., 2014; Brassicales: Barker etal., 2009; Schranz etal., 2011;
Rabosky&Adams,2012;Raboskyetal.,2012)havearguedthat Brassicaceae:Haudryetal.,2013;Kagaleetal.,2014;Asteraceae:
theabsenceof apositive relationship betweenclade age and spe- Barker etal., 2008), but tests of the effect of ancient WGDs on
cies richness calls into question the validity of using birth–death ratesofdiversificationhavenotbeenperformedbeforenow.Sch-
models (and their variants). Drawing on results from simulation ranzetal.(2012)presentedaverbalmodel–theWGDradiation
(Rabosky,2009a,2010),theseauthorsconcludethatheterogene- lag-time model – that hypothesizes that WGDs often result in
ity in diversification rates cannot alone explain the lack of rela- diversificationrateincreases,butfollowingadelay(ofpotentially
tionship between age and richness in empirical data (Rabosky, millionsofyears),resultingintreeimbalancefollowingthedupli-
2009a, 2010; Rabosky & Adams, 2012; Rabosky etal., 2012), cationevent.ThissuggeststhatWGDspromote,butarenotsuf-
and that such a pattern suggests that speciation and extinction ficienttocause,increaseddiversification.
rates are influenced by ecological limits to diversity. These We investigated the link between upticks in rates of diversifi-
authorsgoontoclaimthatmethodsbasedonbirth–deathmod- cation and nine well-documented ancient WGDs that we were
els(e.g.MEDUSA)areinappropriatewhereapositiveage–diver- able toplace onour tree,incorporating uncertainty in theplace-
sity relationship is not supported (Rabosky & Adams, 2012; ment of several of these (Table1), and demonstrated significant
Raboskyetal.,2012). statistical support for a nonrandom association between WGD
However, in a recent simulation study, Stadler etal. (2014) events and a delayed increase in rates of diversification (Fig.3).
examined this question in further detail, demonstrating that Explanations for the lag in diversification rate increases follow-
when stem ages are used to estimate the relationship with diver- ing WGD remain unclear, and could simply reflect early extinc-
sity (as we did in this study), both negative and positive slopes tion events – radiations that do not involve WGD are often
can be produced by a variety of processes, including rate hetero- asymmetrical too (Mooers & Heard, 1997). Nevertheless, Sch-
geneous models. Additionally, Stadler etal. (2014) found that ranz etal. (2012) provide a potential scenario to explain the
the results depend strongly on how higher taxa are delimited, observed lag pattern between WGDs and subsequent upticks in
(cid:1)2015TheAuthors NewPhytologist(2015)207:454–467
NewPhytologist(cid:1)2015NewPhytologistTrust www.newphytologist.com