PDF(915K)

14 Pages·2012·0.89 MB·English

ORIGINAL ARTICLE
doi:10.1111/j.1558-5646.2012.01723.x
INTEGRATING FOSSILS WITH MOLECULAR
PHYLOGENIES IMPROVES INFERENCE
OF TRAIT EVOLUTION
GrahamJ.Slater,1,2 LukeJ.Harmon,3,4 andMichaelE.Alfaro1
1DepartmentofEcologyandEvolutionaryBiology,UniversityofCalifornia,610CharlesEYoungDriveSouth,LosAngeles,
California90095–1606
2E-mail:[email protected]
3DepartmentofBiologicalSciences,UniversityofIdaho,Moscow,Idaho83844
4InstituteforBioinformaticsandEvolutionaryStudies(IBEST),UniversityofIdaho,Moscow,Idaho83844
ReceivedFebruary1,2012
AcceptedMay24,2012
DataArchived:Dryadrepository:doi:10.5061/dryad.q96d7
Comparativebiologistsoftenattempttodrawinferencesabouttempoandmodeinevolutionbycomparingthefitofevolutionary
modelstophylogeneticcomparativedataconsistingofamolecularphylogenywithbranchlengthsandtraitmeasurementsfrom
extant taxa. These kinds of approaches ignore historical evidence for evolutionary pattern and process contained in the fossil
record.Inthisarticle,weshowthroughsimulationthatincorporationoffossilinformationdramaticallyimprovesourabilityto
distinguishamongmodelsofquantitativetraitevolutionusingcomparativedata.WefurthersuggestanovelBayesianapproach
thatallowsfossilinformationtobeintegratedevenwhenexplicitphylogenetichypothesesarelackingforextinctrepresentatives
ofextantclades.Byapplyingthisapproachtoacomparativedatasetcomprisingbodysizesforcaniformcarnivorans,weshow
thatincorporationoffossilinformationnotonlyimprovesancestralstateestimatesrelativetothosederivedfromextanttaxa
alone, but also results in preference of a model of evolution with trend toward large body size over alternative models such
asBrownianmotionorOrnstein–Uhlenbeckprocesses.Ourapproachhighlightstheimportanceofconsideringfossilinformation
whenmakingmacroevolutionaryinference,andprovidesawaytointegratethekindofsparsefossilinformationthatisavailable
tomostevolutionarybiologists.
KEY WORDS: Ancestralstates,Brownianmotion,Caniformia,modelselection,phylogeneticcomparativemethods.
“Iamtryingtopursueasciencethatisbeginningtohaveagood tionsofmacroevolutionaryprocesses.Thestructureofphyloge-
many practitioners but that has no name: the science of four- nies when combined with models of continuous trait evolution
dimensional biology or of time and life. Fossils are pertinent to now allows biologists to address hypotheses about the manner
thisfieldwhentheyaretreatedashistoricalrecords...”(Simpson in which standing phenotypic diversity was generated by test-
1953) ing for shifts in rates of trait evolution (O’Meara et al. 2006;
The last decade has witnessed a revolution in evolutionary Thomas et al. 2006; Eastman et al. 2011; Revell et al. 2012)
biology. Fueled by the increased availability of time-calibrated or adaptive trait optima (Butler and King 2004; Beaulieu et al.
molecularphylogenies(DrummondandRambaut2007)andnew 2012). Trait-dependent diversification models can also be used
analytical tools for analyzing trait data, comparative biologists todeterminetheimpactofthosetraitsonlineagediversification
have increasingly switched their attention from describing evo- (Maddisonetal.2007;FitzJohnetal.2009;FitzJohn2010).Tests
lutionarypatternsandcorrelationstotestingfundamentalpredic- for early bursts (EBs) of phenotypic and lineage diversification
(cid:2)C 2012TheAuthor(s).Evolution(cid:2)C 2012TheSocietyfortheStudyofEvolution.
3931 Evolution66-12:3931–3944
G. J. SLATER ET AL.
consistentwithclassicadaptiveradiationtheory(Simpson1944, nodepriorswhenperformingmacroevolutionaryanalysesonex-
1953; Schluter 2000; Harmon et al. 2003; Gavrilets and Losos tanttaxa.Weshowthatthisapproachperformswell,evenwhen
2009)havebeenaparticularlyrichareaofresearch(e.g.,Pybus prior distributions are not exactly centered on the true values,
andHarvey2000;RaboskyandLovette2008a,b;Harmonetal. as might be the case when sampling from the incomplete fossil
2010;BurbrinkandPyron2010;Mahleretal.2010;Slateretal. record. We finally apply this new approach to the case of body
2010). sizeevolutionincaniformcarnivores(FinarelliandFlynn2006),
Ashistoricalrecordsoftheevolutionaryprocess,fossiltaxa andshowusingamodestnumberofinformativenodepriorsthat
potentiallyprovideatremendoussourceofinformationregarding incorporation of fossil information using our approach not only
thetempoandmodeoflineagediversificationandtraitevolution. improvesparameterestimates,butalsoresultsintheselectionof
However,thecombiningofpaleontologicalandneontological(es- a model with an evolutionary trend toward large body size that
peciallymolecular)dataintophylogeneticcomparativemethods otherwise could not be detected with data derived from extant
hasbeenslowandthefullimpactoffossilinformationonanalyses taxa.
ofextanttaxaremainsrelativelyunexplored.Afewstudieshave
notedexplicitconflictbetweeninferencesaboutlineagediversi-
Methods
ficationratesderivedfrommolecularphylogeniesandthosede-
rivedfromthefossilrecord(QuentalandMarshall2009;Quental MODELS OF TRAIT EVOLUTION
andMarshall2010;Simpsonetal.2011;Rosenblumetal.2012; Inthisarticle,weconsiderfourcommonmodelsthathavebeen
seealsoLiowetal.2010forasimulation-baseddemonstration). appliedtoquantitativetraitdata.Thefirstandmostbasicevolu-
Becausepaleobiologistsandneontologistsusedifferentkindsof tionarymodelforcontinuoustraitsusedbycomparativebiologists
data to infer different aspects of net diversification rates (Foote isBrownianmotion(BM),whichreasonablyapproximatesevo-
1996;Wagner2000;Ricklefs2007;QuentalandMarshall2010), lutionunderdriftoradaptiveevolutiontrackingarandomlyfluc-
this discordance is, perhaps, unsurprising and has primarily led tuatingoptimum(Cavalli-SforzaandEdwards1967;Felsenstein
toeffortstodevelopmathematicalmodelsthatcanbetteraccount 1973, 1985; Hansen and Martins 1996). Under BM, the vector
forunobservedoriginationsandextinctionsinmolecularphylo- of trait values, X, for a clade of n species follow a multivariate
geniesofextanttaxa(e.g.,EtienneandApol2009;Rabosky2009; normaldistribution,
Morlonetal.2011;Etienneetal.2012).
Formodelsoftraitevolution,dataintegrationistheoretically X(1,2,...n) ∼ N(a,σ2C),
morestraightforward.Iffossilspecieswithmeasuredtraitscanbe
placedinaphylogeny,modelsoftraitevolutioncaneasilybefitted whereaisann×1vectorinwhicheachelementistheexpected
to data comprising fossil and extant taxa. It has long been rec- valueofthetrait(which,underBM,istherootstatea),σ2 isthe
ognizedthatincorporatinginformationfromfossiltaxainphylo- Brownian rate parameter and C is an n × n covariance matrix,
geneticreconstructionbasedonmorphologicaldata(Felsenstein referredtoasthephylogeneticvariance–covariancematrix.The
1978;Donoghueetal.1989;Gauthieretal.1988;Huelsenbeck off-diagonalelementsofC,C ,representthesharedpathlength
ij
1991)andinancestralstatereconstructionsofquantitativetraits fromtherootofthephylogenetictreetothemostrecentcommon
(OakleyandCunningham2000;FinarelliandFlynn2006;Alberts ancestor of the ith and jth taxa, whereas diagonal elements C
ii
etal.2009)canimproveparameterestimatesandoverallmethod representthedistancefromtherootnodetotheithterminaltaxon,
performance. The same should be true when comparing the fit whichmaybeextantorextinct.Mostalternativemodelsoftrait
ofdifferentmodelsofquantitativetraitevolutiontocomparative evolution,includingthethreethatweconsiderinthisarticle,are
data, but this has not yet been thoroughly explored. Adapting variants on this basic BM model that involve a transformation
current methods would require knowing both where the fossil of a, σ2C (also referred to as V), or both. Under BM with a
connectstothetreeandthebranchlengthleadingfromthetreeto directionaltrend(Trend:Hunt2006),variancesandcovariances
thefossil.Asignificantbarriertoincorporatingfossilinformation remainthesamebuttheexpectedvalueofthetraitchangesthrough
inanalysesofphenotypicevolution,however,isthatweoftenonly timeasfunctionofatrendparameter,M,thatcanbepositiveor
have very general ideas about their relationships to extant taxa. negative. Under a Brownian process with time-dependent rates,
Inthesecases,integratingfossilinformationrequiresadifferent alsoreferredtoasAC/DC(Accelerating/Deceleratingevolution:
approach. Blombergetal.2003),therateofphenotypicevolutionincreases
Inthisarticle,weexplore,viasimulation,theimpactoffossil (AC)ordecreases(DC)exponentiallythroughtimeasafunction
taxaonmodelselectionwhenfittingmodelsofquantitativetrait of a parameter r that can vary between –∞ and ∞. The DC
evolutiontocomparativedata.WethenpresentanovelBayesian part of this model has also been called EB by Harmon et al.
approachthatallowsfossildatatobeusedtodefineinformative (2010).Finally,underanOrnstein–Uhlenbeck(OU)process,the
3932 EVOLUTION DECEMBER2012
FOSSILS, PHYLOGENIES, AND MODELS OF TRAIT EVOLUTION
Table 1. Parametervaluesandsampledrangesusedforsimula-
tiontests.UnderTrend,theparameterismu,thetrendparameter.
Under SSP, the parameter is α, the strength of selection. Under
AC/DC,theparameterisr,theexponentialchangeparameter.
Model A σ2 non-BMparameter
BM 0 U[0,0.5] –
Trend 0 0.1 U[00.2]
SSP 0 0.1 U[0.01,10]
AC 0 0.0013 U[0,0.1]
DC 0 0.13 U[–0.1,0]
60 50 40 30 20 10 0
Time before present
fossils improve model selection more when extinction rates are
Figure 1. Phylogenetic tree used for simulations. The tree was higher(seeAppendixS1forfulldetails).Wenextsimulatedtrait
generatedunderabirth–deathprocesswithbirthrate=0.1and evolutiononthecomplete377-taxontreeunderthefourmodels
deathrate=0.09.Thereare100extanttaxaand277fossiltaxa. of trait evolution: BM, Trend, AC/DC, and SSP. To explore the
effectsofvaryingmodelparametersonmodelselectionwithor
elements of V are transformed according to their distance from withoutfossilinformation,wesimulated500datasetsundereach
thebasalsplitinthephylogenyandbyaso-called“rubber-band” model, with model parameters drawn at random from uniform
parameter,α(Hansen1997;ButlerandKing2004).UnderOU, ranges (Table 1). The root state was fixed at a = 0 for simula-
as a trait evolves away from its optimal value, it is pulled back tionsunderallmodels.ForBM,wevariedtherateparameterσ2.
towardtheoptimumwithastrengthcorrespondingtoα.Although Forallothermodels,wefixedσ2 andvariedthemodel-specific
itispossiblefortherootstateandtraitoptimumtodiffer,orfor parameterasspecifiedinTable1.Tofullyexploremodel-fitting
multipletraitoptima(ButlerandKing2004)orαvalues(Beaulieu performance,wetreatedAC/DCastwoseparatemodels(ACand
etal.2012)toexist,comparativemethodsoftenassumethatthe DC)forsimulationpurposes.
root state and optimum are the same, reducing OU to a single Foreachsetofsimulations,weinitiallyfitthefourmodels
stationary peak (SSP) model (Harmon et al. 2010) that we will to datasets comprising extant taxa only using maximum likeli-
consider here. Useful summaries of these models and how they hood(ML).Wethenrepeatedthisprocedurefivetimes,eachtime
relate to BM can be found in Hunt (2006) and Harmon et al. randomly adding 5%, 10%, 25%, 50%, and 100% of the total
(2010). numberoffossiltaxatothedataset.Modelfittingwasdoneus-
ing the fitContinuous function in geiger (Harmon et al. 2008).
THE EFFECT OF FOSSIL TAXA ON MODEL SELECTION For each level of fossil sampling, we compared support for the
Wefirstassessedtheimpactofincludingfossiltaxawithknown fourmodelsusingsmall-samplecorrectedAkaikeweights(AIC ;
c
phylogenetic placement and branch lengths on model selection BurnhamandAnderson2002),wherethesamplesizeequalsthe
performance. We generated a phylogenetic tree containing 100 numberofsampledterminalsinthephylogeny,includingfossils
extant taxa under a time-homogeneous, bifurcating birth–death whenapplicable.Iffossilshaveapositiveeffectonmodelselec-
processusingthegeigerpackage(Harmonetal.2008)forR(R tion,thenAkaikeweightsforthetruemodelshouldincreaseas
DevelopmentCoreTeam2011).Atreeofthissizeislargeenough morefossiltaxaareadded.
to potentially allow for differentiation of models of trait evolu- Theadditionoffossiltaxatoadatasetcomprisingextanttaxa,
tion (FitzJohn et al. 2009; FitzJohn 2010; Harmon et al. 2010; asdescribedabove,increasesthetotalnumberoftaxa.Thus,any
Boettiger et al. 2012) and is comparable to one that might be improvementinmodel-fittingperformancecouldbeattributable
used in a typical empirical analysis of extant taxa. To generate toalargernetsamplesize,ratherthanjusttheinclusionoffossil
atreewithalargenumberofextincttips,wesimulatedundera taxa.Wethereforealsoinvestigatedwhether,givenafixedsam-
birth–death process with high relative extinction rate (λ = 0.1, plesizeoftaxa,weshouldprefertosamplefossilsorextanttaxa.
μ=0.09).Theresultingtreehadatotalof377taxa,witharoot Werepeatedthemodel-fittingproceduredescribedabovebutran-
ageof60.95timeunits(Fig.1).Toavoidconfoundingvariation domlysubstitutedproportionsofextanttaxaforfossiltaxa,rather
duetofossilinclusionwithvariationassociatedwithphylogenetic than adding them, thus maintaining a total of 100 taxa. We re-
treeshape,wefocushereonanalysesfromasinglesimulatedtree. peated this procedure for proportions of 5%, 10%, 25%, 50%,
Resultsfromfittingmodelstoalternativephylogenetictreeswere and95%ofthetotalnumberofextantspecies,using95%rather
always qualitatively consistent, although we note that overall, than100%inthefinalstepsothattheresultingtreemaintained
EVOLUTION DECEMBER2012 3933
G. J. SLATER ET AL.
the same root-tip distance without rescaling of branch lengths. Wesubsequentlyperformedrunswhereinformativenormalprior
ModelsupportwasagainassessedusingAIC weights. distributionswereappliedtorandomsamplesof5%,10%,25%,
c
and50%ofnodesinthetree,excludingtherootnode.Foreach
USING FOSSILS AS NODE PRIORS node,themeanofthepriorwassettothetruenodevaluetaken
Formostcomparativebiologists,phylogeneticinformationistyp- fromthesimulateddata,andthestandarddeviation(SD)wasset
ically lacking for fossil members of their study clade. In these atanarbitrarilysmallvalueof0.01.InAppendixS1,wepresent
cases,integratingavailablefossilinformationrequiresadifferent resultsfromanadditionaltestofourapproachusingnodepriors
approachtoML.Underthemultivariatenormaldistribution,the sampledfromasimulatedfossil.Inthesetests,informativenode
log-likelihoodforasetofobservedtraitvaluesXiscomputedas priorsweredefinedbasedonfossiltaxathat(1)weredescendents
⎛ (cid:5) (cid:6)⎞ oftheimmediateparentnodetothenodeofinterestand(2)went
1
⎜exp [X−E(X)](cid:5)(V)−1[X−E(X)] ⎟ extinctwithina10millionyearwindowcenteredonthatnode.As
ln(L)=ln⎜⎝ 2 (cid:7) ⎟⎠, aresultofthissamplingscheme,notallnodespossessedinforma-
(2π)n·det(V)
tivepriors,andthosethatdidhadpriorsthatwerenotnecessarily
centered exactly on their true values. Furthermore, some priors
where E(X) is a n × 1 vector of expected trait values, prime were defined based on only a few fossils with distant phyloge-
specifiesthatthetransposeistaken,andVisthemodel-specific neticrelationshipstothefocalnode(seeAppendixS1formore
variance–covariancematrix(O’Mearaetal.2006;Harmonetal. details).Forbothsetsofanalyses,nodesweresampledatrandom
2010).Normally,onlytraitvaluesfortheterminaltaxaarerepre- amongsimulationsandsamplingschemes.Wesetuninformative
sentedinthisexpression;likelihoodapproachestofittingmodels priorsofU[–∞,∞]ona,ln(σ2)andtheremainingnodesforeach
of trait evolution avoid explicit reconstruction of trait values at model,andonMfortheTrendmodel.ForSSP,weusedauniform
nodesotherthantherootbyintegratingparameterestimatesover prioronαof[10−5,∞].ForAC/DC,weappliedauniformprior
all possible ancestral states (Pagel 1994, 1997; Schluter et al. U[–0.5,0.5]totheexponentialchangeparameter,r.Eachanalysis
1997).InaBayesianframework,thisisequivalenttoassuminga wasrunfor200,000generations,samplingfromthechainevery
uniformpriordensityoneachancestralstate,withtheresultthat 100steps,withthefirst25%ofsamplesdiscardedasburn-in.For
theirposteriordensitiesaredirectlyproportionaltothelikelihood eachfossilsamplingregime,thethreenon-Brownianmodelswere
surface (Maddison 1991; Schluter et al. 1997). However, if we compared to BM using Bayes factors computed from harmonic
lack precise phylogenetic information regarding the placement means of the sampled likelihoods (Newton and Raftery 1994;
fossil taxa butthey canbe associatedwith particular nodesin a Kass and Raftery 1995). To assess false-positive rates, that is,
molecular phylogeny, then their traits could be used to inform casesinwhichweselectanincorrectmodelwithstrongsupport,
priordistributionsplacedonthosenodeswhenfittingmodelsto weassumedthatavaluefor2∗Ln(Bayesfactor)of2orgreater
datasetscomprisingextanttaxa. indicatedpositivesupportforthatmodeloverthenull,BMmodel
WeimplementedaMarkovchainMonteCarlo(MCMC)al- (KassandRaftery1995;Raftery1996).
gorithmtosampletraitevolutionarymodelparametersandnode
valuesfromtheirposteriordistributionswhileallowinginforma-
tivepriordistributionstobeplacedonnodeswithassociatedfossil
Results
taxa.AteachstepinourMCMCsampler,weproposemodelpa-
rametersandancestralstates,andassesstheposteriorprobability THE EFFECT OF FOSSIL TAXA ON MODEL SELECTION
ofthemodelastheproductofthelikelihoodandpriorprobabili- Adding fossil taxa with known phylogenetic placement and
tiesofproposedparametersandancestralstates.Thedecisionto branch lengths had a great impact on model selection perfor-
acceptorrejectproposedparametersaremadeusingthestandard mance, but only when the true model deviated from pure BM.
Metropolis–Hastingsalgorithm(Metropolisetal.1953;Hastings When BM was the true model, we found low power to favor it
1970).Rcodetoperformtheseanalysesiscontainedinthefunc- over other candidate models based on extant taxa only. Adding
tionfitContinuousMCMC(·),whichwillbeavailableviaCRAN fossils had no effect on model selection for any value of σ2
inaforthcomingreleaseofthegeigerlibrary,andisavailablein (Figs. 2A, S6A; Table S3). For the other models, adding fossil
theshorttermasapackagedownloadfromLJH’ssoftwarepage taxaincreasedboththeoverallpowertodetectthetruemodeland
(http://www.webpages.uidaho.edu/∼lukeh/software/index.html). thepowertodosoatsmallparametervalues(i.e.,thoseresulting
Totestthisapproach,weperformedfiveMCMCrunswith in a model most similar to BM) compared to fits using extant
eachsetofsimulations(seeabove)undereachofthefourmodels taxaonly.Aparticularlystrikingresultwasrecoveredwhenthe
ofevolution.Inthefirstsetofruns,thefourmodelswerefitted true model was accelerating evolution. Here, AC was never fa-
tothesimulateddatawithuninformativepriorsonallnodestates. voredoverothermodelsusingextanttaxaonly(Figs.2D,S6M;
3934 EVOLUTION DECEMBER2012
FOSSILS, PHYLOGENIES, AND MODELS OF TRAIT EVOLUTION
ght1.0 A 1.0 B 1.0 C 1.0 D 1.0 E
ei
W5 5 5 5 5
e 0. 0. 0. 0. 0.
k
ai
k0 0 0 0 0
A0. 0. 0. 0. 0.
0.0 0.5 1.0 0.0 0.5 1.0 0.0 0.5 1.0 0.0 0.5 1.0 0.0 0.5 1.0
proportion of available fossils sampled
ght 1.0 F 1.0 G 1.0 H 1.0 I 1.0 J
ei
W 5 5 5 5 5
e 0. 0. 0. 0. 0.
k
ai
k 0 0 0 0 0
A 0. 0. 0. 0. 0.
0.0 0.5 1.0 0.0 0.5 1.0 0.0 0.5 1.0 0.0 0.5 1.0 0.0 0.5 1.0
proportion of 100 taxa that are fossils
BM Trend ACDC SSP
Figure 2. Median Akaike weights for the four evolutionary models fitted by maximum likelihood to 500 datasets. Datasets were
simulatedunderBM(A,F),Trend(B,G),DC(C,H),AC(D,I)andSSP(E,J).Inthetoprow(A–E),increasingproportionsofthe277available
fossils were randomly added to the 100-tip ultrametric tree. In the bottom row (F–J), proportions of extant taxa were replaced with
fossils.
TableS3).Instead,BMwasfavoredatlowincreasesinrateand USING FOSSILS AS NODE PRIORS
SSPathigherincreases(Fig.S6M).Theadditionoffossilsdra- WhenBMwasthetruemodel,therewasverylittlepowertodis-
maticallyincreasedsupportforthecorrectmodelofAC/DCwhile tinguish it from the three alternative models using Bayes factor
simultaneously reducing support for SSP and BM (Figs. 2D, comparisonwithuninformativepriorsonallnodestates.Aswith
S6N–P). MLestimates,BMcouldnotberejectedhowever(Fig.3A).There
Wefoundthatswappingextanttaxaforfossilswhilemain- wasstillpowertodistinguishthethreeothermodelsfromBMus-
tainingaconstanttotalnumberoftaxaalsoresultedinimproved inguninformativenodepriors,particularlywhenparametervalues
model selection performance, relative to analyses based on ex- werelarge(Fig.3G,J,M).WhenSSPwasthetruemodel,itcould
tant taxa only. This result suggests that on a per-taxon basis, bereadilydistinguishedfromBMatallvaluesofα,althoughAC
fossilscontainmoreinformationregardingthemodeoftraitevo- alsoreceivedpositive,albeitlowersupportoverBM.WhenAC
lutionthandoextanttaxa.WhenBMwasthetruemodel,Akaike wasthetruemodel,wefoundittobeindistinguishablefromSSP,
weights for BM remained slightly higher than those for the al- evenatlargevaluesofr,althoughbothreceivedpositivesupport
ternative models, but still ambiguous under all levels of fossil relativetoBM.
swapping (Figs. 2F, S7A–D; Table S4). For all other models, Addition of informative priors on node values improved
swapping fossils for extant taxa both increased support for the modelselectionperformanceinallcasesexceptwhereBMwas
truemodelrelativetoothermodelsandthepowertodiscernthem thetruemodel(Fig.3B,C).BMremainedindistinguishablefrom
atlowparametervalues(Figs.2H–J,S7E–T;TableS4).Oneno- thealternativemodelsatallvaluesofnodesampling,andfalse-
tableexceptionwasfoundwhenthetruemodelwasSSP.Inthis positive rates did not decrease appreciably, except for compar-
case,weightsforSSPdecreasedfromamedianweightof1when isonsofSSPtoBM(Table2).Fortheothermodels,additionof
nofossilsweresampledtoamedianweightof0.59with5%re- priorsonnodevaluesincreasedpowertoidentifythetruemodel,
placementofextanttaxabyfossils(Fig.2J;TableS4).Thisresult withincreasingnumbersofpriorsmarkedlyincreasingpowerat
mayreflecttheobservationsthat:(1)SSPandACaredifficultto lowparametervalues(Fig.3).Themostnotableresultwasrecov-
distinguishwithoutfossilinformation,and/or(2)withoutfossils, eredwhenthetruemodelwasAC.AlthoughACandSSPwere
SSPtendstobefavoredoverAC,regardlessofwhichisthetrue equallyfavoredoverBMbasedonextanttaxaalone,additionof
model. Reducing the number of extant taxa and adding a few fossilinformationintheformofnodepriorsresultedinincreased
fossilsseemstohavetheeffectofaccentuatingthisuncertainty. supportforthetruemodel(AC),whilesupportforSSPdecreased.
However, increasing the proportion of fossils swapped to 10% Similarly,increasingthenumberofnodeswithinformativepriors
increased median weights for SSP to 0.99, suggesting that high increasedsupportforSSPwhiledecreasingsupportforACwhen
levelsoffossilsamplingarenotnecessarytoincreaseconfident SSPwasthetruemodel.InAppendixS1,weshowthatourap-
identificationofthetruemodel. proachalsoperformswellusinginformativepriorsderivedfrom
EVOLUTION DECEMBER2012 3935
G. J. SLATER ET AL.
A B C
3 3 3
F) 2 2 2
B
n( 1 1 1
2 * l 0 0 0
1 1 1
- - -
0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8
2 2 2
60 D 60 E 60 F
F) 0 0 0
B 4 4 4
n(
2 * l 20 20 20
0 0 0
0.02 0.04 0.06 0.08 0.02 0.04 0.06 0.08 0.02 0.04 0.06 0.08
0 0 0
2 G 2 H 2 I
1 1 1
BF) 80 80 80
n(
2 * l 40 40 40
0 0 0
-0.08 -0.04 -0.08 -0.04 -0.08 -0.04
r r r
0 J 0 K 0 L
3 3 3
F)
B 0 0 0
n( 2 2 2
2 * l 10 10 10
0 0 0
0.02 0.04 0.06 0.08 0.02 0.04 0.06 0.08 0.02 0.04 0.06 0.08
r r r
M N O
F) 0 0 0
n(B 15 15 15
2 * l 0 0 0
5 5 5
0 0 0
2 4 6 8 2 4 6 8 2 4 6 8
Figure 3. Median Bayes factors for candidate models fitted to 500 simulated datasets with informative node value priors on 0, 10,
and50%oftotalnodesinthetree.Ineachcase,candidatemodelsarecomparedtoanullmodelofBrownianmotion.Rowsrepresent
simulations under BM (A–C), Trend (D–F), DC (G–I), AC (J–L), and SSP (M–O). In each panel, the x-axis represents the range of model
parametervaluessimulatedunder,whereasplottedlinesrepresentedthemediansofbinnedBayesfactors.LinecolorsareTrend(gray),
AC/DC(blue),andSSP(black),allcomparedtoBM.ThedashedredlineindicatestheBayesFactorrequiredtoachievepositiveevidence
infavorofagivencandidatemodel.
asimulatedfossilrecord,despitethefactthatpriorswerederived iformcarnivores(FinarelliandFlynn2006).Morphologicaland
fromasmallnumberoffossilsthatwereoftenyoungerthatthe ecologicaldiversitywithintheCaniformiaisparticularlystriking;
nodeofinterest. caniformsrangefromexclusiveherbivorestospecialistcarnivores
(Figueiridoetal.2010)andspantheentirerangeofbodysizesen-
BODY SIZE EVOLUTION IN CANIFORM compassedbyextantCarnivora(overfourordersofmagnitude).
CARNIVORANS FinarelliandFlynn(2006)demonstratedthatunreasonablylarge
ToillustratetheapplicationofourBayesianapproachtoanem- estimates(∼25kg)fortheancestralbodysizeofcaniformswere
piricaldataset,werevisitedthecaseofbodysizeevolutionincan- obtained when only extant taxa are considered. However, when
3936 EVOLUTION DECEMBER2012
FOSSILS, PHYLOGENIES, AND MODELS OF TRAIT EVOLUTION
Table 2. False-positiveratesforBayesianmodelselectionusingnodevaluepriorscenteredonthetruenodevalue.Falsepositivesunder
threeBayesfactorcut-offsaregiven.
2×lnBF=2 2×lnBF=6 2×lnBF=10
Prop.nodes Trend AC/DC SSP Trend AC/DC SSP Trend AC/DC SSP
0 0.00 0.11 0.11 0.01 0.01 0.01 0 0 0
0.05 0.09 0.11 0.07 0.01 0.01 0.01 0 0 0
0.1 0.09 0.12 0.07 0.01 0.01 0.01 0 0 0
0.25 0.09 0.11 0.07 0.01 0.01 0.01 0 0 0
0.5 0.10 0.12 0.07 0.01 0.01 0.01 0 0 0
fossil taxa were included in their phylogeny, more reasonable
Procyonidae
ancestral size estimates of ∼5 kg were obtained. Finarelli and
Flynn (2006) suggested that this might be the result of an evo-
Mustelidae
lutionary trend toward large body size (Cope’s rule, e.g., Alroy
1998;PurvisandOrme2005)insomeextantfamilies.Wethere-
Ailuridae
fore asked whether the incorporation of fossil information into Mephitidae
acomparativedatasetderivedfromextantcaniformcarnivorans Phocidae
couldprovidesupportforthismodeloverothercandidatemodels
Otariidae
ofbodysizeevolution.
Odobenidae
Ursidae
CANIFORM METHODS
Canidae
Bodymassdataforextantcaniformspecieswereobtainedfrom
thepanTHERIAdatabase(Jonesetal.2009),withmissingspecies
addedfromtheliteraturewherepossible.Weinitiallyattemptedto 50 40 30 20 10 0 1 100 2000
Millions of Years Ago body mass
usethecarnivoranportionofthemammaliansupertree(Bininda-
(kg)
Emondsetal.2007).However,attemptstoidentifyfossiltaxathat
couldbephylogeneticallyandtemporallyassociatedwithnodes Figure 4. Time-calibratedphylogenyof135speciesofextantcan-
iformsandassociatedbodymasses(notethatbarlengthsarerep-
in that tree revealed some significant discrepancies in topology
resented on a log scale). Colors indicate caniform families. Gray
and divergence times compared with current understanding of
circlesatnodesindicatelocationsoffossilpriorsusedforfitting
carnivoranphylogeny(Bardelebenetal.2005;Flynnetal.2005;
modelsofbodysizeevolution.Theblacknodelabelindicatesthe
Koepflietal.2007;Koepflietal.2008;Krauseetal.2008;Eizirik rootprior.
etal.2010;FultonandStrobeck2010a,b;Meredithetal.2011).
We therefore constructed a new time-calibrated phylogeny for
caniforms by grafting published and new time-calibrated phy- The final set of priors spanned 11 internal nodes, plus the root
logenies for caniform clades onto a family-level backbone tree node(Fig.4;Table3).Bodymassesforindividualfossiltaxaused
(Eiziriketal.2010).Fulldetailsofphylogenyreconstructionare todefineinformativenodepriorsareprovidedinTableS7.
providedinAppendixS2and.xmlfileshavebeendepositedon Wefirst comparedthefit ofthree modelsoftrait evolution
Dryad.Afterremovingnonoverlappingspeciesfromthetreeand (BM,AC/DC,SSP)withuninformativepriorsplacedonallnode
bodysizedata,weretainedafinalsampleof135taxa,approxi- values.Twochainswererunfor500,000generationseach,sam-
mately82%ofextantcaniformdiversity(Fig.4).Wesubsequently plingfromthechainevery100steps,andthefirst10%ofsamples
surveyedtheliteratureforfossiltaxathatcouldbeassociatedwith fromeachrunasburn-inbasedoninspectionoflikelihoodtrace
nodes in the caniform phylogeny, where possible giving prefer- plots.Marginallikelihoodsandparametervalueswerecomputed
encetofossilsthathadbeensubjecttophylogeneticanalyses.We from the retained sample. We then re-ran analyses using the
usedpublishedestimatesofbodymasswhenavailable,butother- same MCMC procedures but with informative normal prior
wiseestimatedbodymassesforfossiltaxafromlowerfirstmolar distributionsplacedonthe11internalnodeswithassociatedfossil
lengths,usingregressionequationsprovidedinVanValkenburgh information.Fortheseanalyses,theTrendmodelwasaddedtothe
(1990).Informativenodevaluepriorswerederivedfrommeans suite of models fitted. Because fossil information was available
and SDs computed from natural log-transformed fossil masses. for the root node of caniform phylogeny, we performed a third
EVOLUTION DECEMBER2012 3937
G. J. SLATER ET AL.
Table 3. Priorsusedforfittingmodelsoftraitevolutiontothe HPD:–0.002–0.063),whichisapproximatelyequivalenttoanet
caniformdataset.Meansandstandarddeviationsareinunitsof expectedincreaseinsizeof∼7kgoverthe50millionyearhis-
naturallog(Ln)kilogram. tory of extant caniforms. Based on simulations, a trend of this
magnitude is difficult to favor with positive evidence or greater
Node Mean SD
withinformativepriorsononly10%ofnodes(Fig.3E;Appendix
Rootnode 0.749 0.485 S1).ThefactthatwefindpositiveevidenceinfavorofTrendwith
MRCACaninae 1.451 0.487 only8%ofavailablenodessampledisthereforeencouraging.It
MRCAVulpini 1.974 0.29
is also worth noting that the Trend model used here assumes a
MRCAVulpes 1.831 0.211
homogenous directional trend across all lineages through time
MRCACanini 2.225 0.184
(Hunt2006).Thismightbeexpectedunderaphenomenonsuch
MRCAArctoidea 0.615 0.336
as Cope’s rule where an active trend drives all clade members
MRCAUrsinae/Tremarctinae 4.293 0.129
MRCAUrsinae 4.406 0.028 toward larger size (Alroy 1998; Van Valkenburgh et al. 2004;
MRCAcrownMustelidae 1.296 0.4645 Purvis and Orme 2005). However, it is possible that the evolu-
MRCAcrownProcyonidae 0.417 0.172 tionofcaniformbodysizesisbetterexplainedbyapassivetrend,
MRCAPinnipedimorpha+Musteloidea 1.14 0.33 suchasdiffusionawayfromareflectiveorabsorbinglowerbound
MRCAcrownPinnipedia 4.111 0.957
(Stanley 1973; Gould 1988; Alroy 2000). Although such an ar-
gumentmakessenseonphysiologicalgrounds(Tracy1977;West
etal.2002),itcannotcurrentlybetestedusingphylogeneticcom-
setofanalyseswithpriorsonthe11internalnodesplustheroot parativemodels.
node.
Discussion
CANIFORM RESULTS
Fossilshadasignificantimpactonrootstateestimation(Fig.5). Fittingmodelsofquantitativetraitevolutiontophylogeneticcom-
With uninformative priors on all internal nodes, we estimated parativedataprovideapowerfulwayoftestingfundamentalhy-
a large body mass (∼18.5 kg) for the ancestral caniform under potheses about macroevolutionary pattern and process (Hansen
allmodels,albeitwithwide95%highestposteriordensityinter- 1997; Blomberg et al. 2003; Butler and King 2004; O’Meara
vals (HPD: 1.85–171 kg; Fig. 5). Fitting the same models with et al. 2006; Thomas et al. 2006; Thomas et al. 2009; FitzJohn
informative prior distributions on 11 internal nodes resulted in 2010;Harmonetal.2010;Eastmanetal.2011;Revelletal.2012;
substantially smaller modal estimates (∼2.5 kg) and narrower Slateretal.2012).Alltoofrequentlyhowever,analysesarelim-
95% HPD intervals for the root state (0.34–10 kg; Fig. 5). Ad- itedtoextanttaxa,neglectingthehistoricalinformationthatfossil
dition of a root node prior unsurprisingly resulted in narrower datacanprovide(Simpson1944,1953;Losos2011).Theresults
still95%HPDintervals(0.89–4.5kg),althoughthemodedidnot presentedinthisstudysuggestthatinclusionofjustafewfossil
shift.Overall,theseresultsareconsistentwiththoseofFinarelli taxanotonlyincreasesthestatisticalpoweroftheseapproaches,
and Flynn (2006) and appear to suggest a non-Brownian model butcanalsoleadtodifferentiationofmodelsthatsimplycannot
ofbodysizeevolutionincaniformcarnivorans. be distinguished based on extant taxa data alone. Furthermore,
Fossilsalsohadanimpactonmodelselectionforthecani- althoughmanycladeslackawell-sampledfossilrecord(apoint
form body size dataset. With uninformative priors on all nodes, wewillreturntoattheendofthisarticle),wehaveshownhow
BM,AC/DC,andSSPwerebroadlyindistinguishable(Table4). sparsefossilinformationwithoutphylogeneticresolutioncanbe
Theadditionoffossil-derivedpriorson11internalnodeshadlittle integratedintomacroevolutionaryanalysisinawaythatcanstill
effect on comparisons between these three models, but did lead improvemodelselectionoverestimatesbasedonextanttaxaonly.
topositive(2×lnBF>2)supportforTrendoverBM.Addition
of a prior on the root state had little impact on model selection DO FOSSILS IMPROVE ESTIMATES OF TRAIT
compared with results obtained using internal node priors only EVOLUTION?
(Table4). We found that fossils had little effect on model selection when
Although we find positive evidence for a trend toward in- BM was the generating model of evolution. Instead, the great-
creasingbodysizeincaniformevolution,supportforthismodel est impact of fossils is realized when the generating model is a
over BM is relatively low. One possible reason for this is the time-heterogeneousprocesses.Onepossibleexplanationforthis
interactingeffectofthestrengthofthetrendandtheproportion resultisthatourmodelcomparisonmetric,Bayesfactorsderived
of nodes for which we had available prior information. Our es- fromtheharmonicmeanestimator(HME)ofthemarginallikeli-
timates for the Trend parameter were low (M = 0.0296, 95% hood,isapoormodelselectiontool.Recentworkindemographic
3938 EVOLUTION DECEMBER2012
FOSSILS, PHYLOGENIES, AND MODELS OF TRAIT EVOLUTION
A B
sity 0.6 sity 0.6
n n
e e
D D
0 0
0. 0.
-2 0 2 4 6 8 -2 0 2 4 6 8
Ln(Mass) Ln(Mass)
C D extant
fossils
sity 0.6 sity 0.6 fossils + root
n n
e e
D D
0 0
0. 0.
-2 0 2 4 6 8 -2 0 2 4 6 8
Ln(Mass) Ln(Mass)
Figure 5. PosteriordistributionsfortherootnodestateinthecaniformexampleestimatedunderBM(A),Trend(B),AC/DC(C),andSSP
(D).ColorsrepresentMCMCanalyseswithandwithoutfossilinformationincorporatedasinformativenodepriors.
Table 4. Resultsofmodelfittingforthecaniformdatasetwithnopriors,withnodepriors,andwithnodeandrootpriors.Marginal
likelihoodsaregivenforallmodelsexceptTrendwithnopriors,whichcannotbedistinguishedfromBM.Bayesfactorsarecomputedfor
modelsrelativetoBM.AsterixindicatespositiveevidenceinfavorofthecandidatemodeloverBM.
Nopriors Nodepriors,norootprior Nodepriors,rootprior
MarginalLk 2×LnBF MarginalLk 2×LnBF MarginalLk 2×LnBF
BM −169.61 – −178.72 – −178.36 –
Trend – – −177.33 2.77∗ −177.02 2.69∗
SSP −169.54 0.13 −178.74 −0.05 −178.64 −0.55
ACDC −169.83 −0.44 −179.17 −0.90 −178.73 −0.74
andphylogeneticmodelselectionhasdemonstratedthatthiscan fossilinformation(itshouldbenotedhoweverthatastrongprior
sometimes be the case (Xie et al. 2011; Baele et al. in press). on the root state would be sufficient to allow this model to be
Alternative metrics, such as AICM (Raftery et al. 2007), and fitted;OakelyandCunningham2000;seealsoFelsenstein1988).
pathsampling(thermodynamicintegration:LartillotandPhillipe Theimprovedpowerfortheothermodelsfittedhereissimilarlya
2006) and stepping-stone (Xie et al. 2011) algorithms can out- functionoftheirtimedependency.Itwouldbeinterestingtofur-
perform HME estimators in many cases (Baele et al. in press). therinvestigatetheimpactoffossilswhenfittingmodelsallowing
The integration of Bayesian approaches into phylogenetic com- forcomponentsofpunctuated(i.e.,speciational)andgradual(i.e.,
parativemethodsisstillinitsinfancy,butfurtherworkonmodel Brownian)evolution(Hunt2007b,2008;Bokma2008).
selection procedures has the potential to greatly increase power Wefoundthatwithoutfossils,thepowertodiscerndecelerat-
todiscernamongmodelsofquantitativetraitevolution(seealso ingevolutionwasgreaterthanthatforacceleratingevolution.This
Boettigeretal.2012). makes sense given the differing expectations of these opposing
The most obvious impact of adding fossil data was for the processes.Deceleratingevolutionresultsinphenotypicdisparity
caseoftheTrendmodel,wheretheexpectedvalueofatraitin- beingpartitionedamongratherthanwithinclades,apatternthat
creasesordecreasesasafunctionoftime.Identificationofatrend canbereadilyidentifiedfromextanttaxadatawhenthedifference
isthereforeexplicitlydependentonknowledgeofboththestart- betweeninitialandendratesislargeenough(Harmonetal.2003).
ing and end states. If a directional trend exists in comparative Adding fossils provides additional information about the tempo
data,ourabilitytoaccuratelyinfertheancestralstateforthetrait atwhichdisparityevolvesandthusfurtherallowsfordifferenti-
fromdataderivedfromextanttaxaonlyislostasthetraitevolves ationofslowdecelerationfromaconstantratesprocess(Figs.2,
away from its initial state (e.g., Oakley and Cunningham 2000; 3, S3, S6, S7). Under accelerating evolution, in contrast, faster
Finarelli and Flynn, 2006; Albert et al. 2009). As a result, we evolutiontowardthetipsresultinincreasedwithin-cladedispar-
completelylacktheabilitytofittheTrendmodelwithoutsome ity relative to among clade disparity—a pattern also seen under
EVOLUTION DECEMBER2012 3939
G. J. SLATER ET AL.
theSSPmodel.Ourresultssuggestthatacceleratingevolutionis usedhere,incontrast,makesexplicituseofsynapomorphy-based
difficulttodistinguishfromBMwithafastratebasedonextant assessmentoffossilplacementthatismorefamiliartoevolution-
taxaaloneunlessthemagnitudeofaccelerationisrelativelylarge. arybiologistsfromdivergencetimeestimation(Drummondand
Insuchcases,however,acceleratingevolutionbecomesdifficult Rambaut2007),andhasmoreincommonwithapproaches,such
todistinguishfromSSP. as ancestor-descendent comparison, that have been traditionally
Fossilinformationdramaticallyimprovedourabilitytodis- beenusedinpaleobiologicalinference(Gingerich1987;Hulbert
tinguishSSPandAC.Acceleratingevolutionresultsinrapidevo- 1993; Wagner 1996; Alroy 1998; Van Valkenburgh et al. 2004;
lutionarychangealongbranchesclosertothetipsofaphylogeny, Honeetal.2005).Ourresultsshowthatmostpowerisobtained
withthenetresultthatrecentspeciespossessmarkedlydifferent when fossils are placed in as precise a phylogenetic context as
traitvaluesthantheirsisterspecies.Lookingatextanttaxaonly, possible. If this is not possible, then using fossils as priors on
theoutcomeofthisprocessisverysimilartoanOUprocess,as node states is still preferable over approaches that ignore fossil
bothtendtoerasephylogeneticsignal(Revelletal.2008).Thein- informationentirely.
corporationofdatafromfossiltaxaprovidesasourceofhistorical A pervasive concern in the use of Bayesian approaches is
information regarding character change toward the root the tree that inappropriate or incorrect priors can dominate theanalysis,
andthusenablesthesetwomodelstobedistinguished.Thisresult leadingtomisleadingsignal.Carefulselectionoffossilstempo-
isparticularlyinterestinggiventhatHarmonetal.(2010)found rallyspanningeithersideofanodeofinterestshouldameliorate
aprevalenceofsignalforSSPoverBMandDC(or“EB”)ina any tendency toward biased or incorrect priors (see Supporting
large number of comparative datasets. Those authors suggested Information).Itshouldalsobenotedthatourapproachallowsfor
thattheEBpatternexpectedunderadaptiveradiationhypotheses the same fossil to be used in defining more than one prior. For
(Simpson1944,1953)mightberareincomparativedataandthat example,afossilmightbetheoldestmemberofacrowngroup
stabilizing selection might be a more pervasive pattern. Adding whilealsobeingastemrepresentativeofamoreinclusiveclade.
fossilinformationisunlikelytoswingapreferredSSPmodelto Inthiscase,thesamefossilcouldbeusedwhendefiningthedistri-
DC; however,weare unabletorule outthat a pervasivepattern butionsplacedonbothnodes.Alogicalextensionofourmethod
of accelerating evolution is obscured in these datasets. Acceler- wouldbetousethetemporalproximityofafossilorsetoffossils
ating evolution might occur in adaptively radiating clades that todeterminethestrengthofthepriorplacedonanode.Usingin-
have yet to saturate ecomorphospace, for example, as a young formationinthiswayshouldleadtomoreappropriatepriorsand
claderadiatestoreplaceawaning,ecologicallyequivalentclade makesmorecompleteuseofpotentiallysparseavailabledata.
(e.g., Van Valkenburgh 1999). If diversifying clades are charac- One limitation of our approach is that fossil diversity from
terizedbythesekindsofpatterns,thenfossilinformationisofan radiationsalongstemlineageswillbemostlyunusableiftheentire
evengreatersignificanceinprovidinginsightintothedynamics stemgrouphasgoneextinct.Forexample,althoughwewereable
ofadaptiveradiation. toassociatefossiltaxawith12nodesinthephylogenyofextant
caniforms, we were unable to include information from the ex-
USING FOSSILS AS NODE PRIORS tinct canidsubfamilies HesperocyoninaeandBorophaginae,the
Wealsoshowedthattraitdatafromfossiltaxacanbemeaning- entire family Amphicyonidae, and the majority of fossil mem-
fullyincorporatedintocomparativeanalysesasinformativeprior bers of the extant, monotypic families Odobenidae (walruses)
distributionsonnodevalues.Thisisparticularlyappealingasfos- andAiluridae(redpanda),despitethefactthatallpossessgood-
silrepresentativesofextantcladesareinfrequentlyplacedinan to-excellentfossilrecords(Deme´re´1994;Wang1994;Hunt1998;
explicit phylogenetic context. Moreover, even when fossils can Wang et al. 1999; Morlo and Peigne´ 2010). Given that biologi-
beplacedinaphylogeny,branchlengthsaretypicallynotavail- callyreasonablevaluesformostevolutionarymodelparameters
able (but see Pyron 2011), rendering fitting of time-dependent arelikelytobeatthelowerendoftheirrespectivescales,wepo-
processesdifficult.Felsenstein(2002)suggestedanapproachby tentiallysurrendertremendousstatisticalpowerbynotincluding
whichfossiltaxacouldbeplacedinamolecularphylogenyand thisavailableinformation.Reassuringly,bothoursimulationsand
theirbranchlengthsinferredusingtheircontinuoustraitvalues.In ourempiricalexampledemonstratethatifatrendintraitvalues
thisscheme,fossilswouldbeplacedbymaximizingthelikelihood or time-dependency of rates actually exists, we retain substan-
of the BM process given knowledge of the phylogeny of extant tial power to detect this with a sample of extant taxa and a set
taxa,asmightbeinferredfrommoleculardata.AlthoughFelsen- ofinformativenodepriorsderivedfromsuitablefossildata.For
stein described this algorithm with a view toward phylogenetic more complex models of trait diversification however, explicit
inferenceratherthanmodelfitting(aformof“TotalEvidence”: phylogenetic hypotheses or alternative approaches for integrat-
Felsenstein 2002, p. 32), it could potentially be applied to the ing phylogenetically uncertain fossil radiations are likely to be
problemthatwehaveattemptedtosolve.Theapproachwehave necessary.
3940 EVOLUTION DECEMBER2012