跳到主要内容
广告
浏览主题
?

点击PLOS分类法找到你所在领域的文章。

有关公共科学图书馆学科领域的更多信息,请单击在这里

  • 加载指标

FastTree 2 -用于大对齐的近似最大似然树

  • 摩根·n·普赖斯

    MorganNPrice@yahoo.com

    从属关系美国加利福尼亚州伯克利市劳伦斯伯克利国家实验室物理生物科学部,美国加利福尼亚州伯克利市劳伦斯伯克利国家实验室微生物压力与生存虚拟研究所

  • 帕拉米尔·s·迪哈尔,

    从属关系美国加利福尼亚州伯克利市劳伦斯伯克利国家实验室物理生物科学部,美国加利福尼亚州伯克利市劳伦斯伯克利国家实验室微生物压力与生存虚拟研究所

  • 亚当·阿金

    从属关系美国加利福尼亚州伯克利劳伦斯伯克利国家实验室物理生物科学部,美国加利福尼亚州伯克利劳伦斯伯克利国家实验室微生物压力与生存虚拟研究所,美国加利福尼亚州伯克利加州大学生物工程系,美国加利福尼亚州伯克利

摘要

背景

我们最近描述了FastTree,这是一种用于推断多达数十万个序列的系统发育比对的工具。在这里,我们描述了FastTree在不牺牲可伸缩性的情况下提高准确性的改进。

方法/主要结果

FastTree 1使用最近邻交换(NNIs)和最小进化准则来改进树,FastTree 2增加了最小进化子树修剪-再嫁接(SPRs)和最大似然NNIs。FastTree 2使用启发式来限制搜索更好的树,并估计每个站点的进化速度(“CAT”近似)。尽管如此,对于模拟和真实的对齐,FastTree 2比最大似然NNIs的标准实现(默认设置的PhyML 3)略微更准确。虽然FastTree 2不像使用最大似然SPRs的方法那样准确,但大多数不一致的分割都得不到很好的支持,对于大型对齐,FastTree 2要快100 - 1000倍。FastTree 2在22小时和5.8 gb内存的台式计算机上推断了237,882个不同的16S核糖体rna的拓扑结构和基于可能性的本地支持值。

结论/意义

FastTree 2允许对巨大的排列进行最大似然系统发育的推断。FastTree 2是免费的http://www.microbesonline.org/fasttree

简介

从相关的DNA或蛋白质序列家族推断进化关系或系统发育是计算生物学的核心方法。基于序列的系统发育被广泛用于了解生物的进化关系和分析基因的功能。

最大的基因家族已经包含了数万到数十万个代表基因,随着DNA测序的快速改进,我们预计更大的数据集很快就会到来。大型家族可以使用基于剖面的方法进行排列,该方法随序列数量线性扩展(http://hmmer.janelia.org/[1]).然而,大多数从这些序列推断系统发育的方法都是O(),或者更糟的是,在哪里序列数是和吗是对齐的长度(宽度)。因此,推断系统发育已经成为计算上的挑战。

我们最近描述了一种可扩展的推断系统发育的方法FastTree 1.0[2].FastTree 1.0基于“最小进化”原则——它试图找到一种拓扑结构,使进化量或分支长度的总和最小化。FastTree 1.0使用了邻居连接的启发式变体[3][4]快速查找起始树,并使用NNIs (nearest-neighbor interchanges)对拓扑进行细化。(最近邻交换交换一个节点和它的邻居;例如,它可能会改变((A、B), C, D), ((A、C), B, D)或((A, D), B, C))。FastTree在O()空间,其中通过存储子树的配置文件而不是它们之间的距离来表示字母表中的字符数。这比存储成对距离所需的内存要少得多,而成对距离对于邻居连接和相关方法是必要的。这也允许将理论运行时间减少到O().(FastTree 1.0也包含了一些O()步骤,但这些已被删除,见http://www.microbesonline.org/fasttree/ChangeLog)。相比之下,计算所有的成对距离,这是大多数最小进化方法所需要的,需要O()时间。与其他最小进化方法相比,FastTree 1.0的主要限制是它不能在初始邻居连接阶段纠正多次替换的距离。然而,NNIs远远弥补了这一点。在实践中,FastTree 1.0比大多数其他最小进化方法更准确,但不如最大似然方法准确[2]

在最大似然(ML)方法中,进化被明确地用转移速率矩阵建模,而最能解释数据的树——具有最高似然的树——就是最好的树[5].ML标准对树进行排序,但没有指定如何找到一个好的拓扑。因为ML系统发育推理是NP完全的[6],没有任何实用的方法可以保证它能找到一个大对齐的最优拓扑。最具可伸缩性的ML方法,如PhyML和RAxML,从由更快的方法生成的起始树开始,并尝试通过优化单个分支长度和执行局部重新排列来增加可能性[7]- - - - - -[9].通过每次移动只重新优化几个分支长度,考虑或执行一个移动的成本可以减少到O()时间,地点是字母的大小。然而,在实践中,移动的数量大致增长为O(),并且优化步骤本身就很慢,因为它们需要数值求解和迭代。这就解释了为什么PhyML和RAxML都需要超过一天的时间来处理1000个蛋白质序列。用自举法估计树的可靠性[10]通常会将计算需求再增加100倍(尽管可以通过跨复制重用计算来减少这一需求)[11]).

在这里,我们将介绍FastTree 2,这是一个用于推断大型对齐ML树的工具。FastTree 2除了使用邻居连接构造初始树并使用最小进化NNIs对其进行改进外,还使用了最小进化子树修剪-再嫁接(SPRs)。[8][12]和ML NNIs来进一步改进树。(在子树修剪-重新嫁接中,子树从树中移除并重新插入到其他地方,例如,修剪和重新嫁接C可能会从((a,B),(C,D),E)变为((a,(B,C)),D,E)。)FastTree 2使用启发式方法来减少搜索空间,从而保持两个阶段的可伸缩性。减少搜索空间的另一个理由是,密集的树搜索通常会发现树的长度或可能性有微小的改进,但这些变化可能在统计学或生物学上不显著(例如,[13]).简单地说,FastTree的关键启发是:

  • 它使用“线性SPRs”来考虑可能的SPR行动。在每个节点上,它首先检查最短的spr,然后扩展最有希望的候选对象。
  • 它为每个子树只搜索两次SPR移动,而不是迭代直到收敛。
  • 在ML阶段,它将ML NNIs限制为最多轮;在实践中,它在达到这个极限之前收敛,但是这个极限确保了可预测的运行时间。
  • 它限制了优化模型参数和分支长度的工作。
  • 它放弃了对NNI移动的优化,在部分优化后,似乎显著降低了可能性。
  • 它并没有试图改进树中最近几轮没有改进的部分。

为了解释不同地点的比率差异,FastTree使用了“CAT”近似[14]而不是标准的四个速率的离散伽马模型([15].有些站点比其他站点进化得慢得多,解释这种情况的理想方法是将每个站点的可能性与该站点的(未知的)相对进化速率相结合,使用先验分布与相对速率相结合,例如gamma分布。然而,这些积分是难以分析和计算的。““方法是使用四个速率类别来近似连续伽玛分布。然而,仍然需要四倍多的CPU时间和内存在不同的站点没有速率变化的模型。此外,对于大型对准,数据严格限制了每个站点的速率。因此,对每个地点的速率(CAT)进行准确的估计,要比对四个潜在速率([14].FastTree从20个固定的可能性中选择每个站点最可能的概率。

由于启发式,FastTree 2不能保证在树空间中达到局部最优似然。然而,在每一步它都保证可能性增加(在CAT近似下)。因此,FastTree 2是一种近似最大似然方法。

我们将在实践中证明,FastTree 2比最大似然NNIs的标准实现(默认设置下的PhyML 3)稍微准确一些[16][17].具体来说,在模拟中,FastTree 2恢复了更高比例的真实分割,在真正的对齐上,FastTree 2的拓扑结构往往有更高的可能性。FastTree的最小进化SPR移动使其比PhyML的起始树更好,PhyML的起始树是用BIONJ(邻居连接的加权变体)获得的[18]).这远远弥补了FastTree的启发式,它降低了ML NNIs的搜索强度,但对准确性几乎没有影响。我们也证实了用CAT近似代替模型(它本身是连续伽玛分布的近似值)对树的质量影响很小。

虽然FastTree 2的精度明显低于使用SPR移动的ML方法,例如设置较慢的PhyML或RAxML,但大多数不一致的分割都不受支持,而FastTree要快得多。FastTree 2可以在台式计算机上在一天内分析数万或数十万个序列的比对。对于500个或更多序列的比对,FastTree 2比PhyML 3.0或RAxML 7.2.1至少快100倍。FastTree 2比RAxML 7快,主要是因为不太密集的ML搜索(NNIs而不是SPRs),而且RAxML 7优化了分支长度模型。但是,FastTree也有一个更快的开始树,它最初比RAxML 7更快地增加可能性。

因为它的速度,FastTree 2适合于引导。然而,为了更快地估计树的可靠性,FastTree 2提供了基于Shimodaira-Hasegawa (SH)测试的局部支持值[16][17][19].FastTree 2对于重建生命树和分析数以百万计的通过基因组测序确定的未特征蛋白质是有用的。

结果

我们将FastTree的速度和准确性与最流行的最大似然方法PhyML 3.0和RAxML 7进行了比较。为了测量结果树的质量,我们测量了模拟对齐的拓扑精度和真实生物对齐的可能性。

仿真中的拓扑精度

我们用250到5000个序列对FastTree进行了模拟蛋白质比对测试[2].这些模拟来自于基因组规模研究中出现的不同基因家族(“同源群集合”或cog,[20]).模拟包括不同地点的不同进化速率,包括间隔的实际位置。模拟可从FastTree网站(http://microbesonline.org/fasttree/#Sims).

我们将拓扑精度定义为每种方法恢复的真实树中分裂的比例。这是拓扑(“Robinson-Foulds”)距离的逆,缩放到0到1的范围。如表1FastTree 2比默认设置的PhyML 3 (NNI搜索)更准确,比最小进化或简约方法更准确,但不如使用SPR移动的ML方法准确。FastTree 2与其他方法的准确性差异有统计学意义(全部,配对测试)。

为了测试使用ML SPR移动发现的额外真分割的实际意义,我们检查了PhyML 3报告的局部支持值。我们将“强支持”定义为同时具有sh样局部支持和近似似然比检验(aLRT)支持[21]95%或更高。在PhyML 3通过SPR移动发现而FastTree 2没有发现的真正分裂中,只有16%得到了强有力的支持。支撑值的完整分布如图所示图1.相反,在PhyML 3用SPRs而不是FastTree发现的强烈支持的分裂中,20%是不正确的。因此,很少有额外的真正分裂有高支持度,而在不同意的分裂中,即使有高支持度的分裂也有很大概率是不正确的。

缩略图
图1所示。本地支持值的分裂发现的PhyML与SPR移动和/或FastTree。

我们检查了由PhyML 3.0推断的分割的局部支持值+ SPRs模拟校准与250个蛋白质序列。我们将PhyML的分裂分类为正确的并被PhyML和FastTree发现,正确的但被FastTree遗漏,或不正确的。我们展示了每个类的支持值的分布。最右边的bin包含了强支持的分割(0.95到1.0),灰色虚线显示了均匀分布。支持值为PhyML近似似然比检验的最小值[21]和SH-like[16][17]当地的支持。

https://doi.org/10.1371/journal.pone.0009490.g001

为了理解为什么FastTree 2在NNI搜索方面优于PhyML 3,我们以FastTree的最小进化树作为它的起始树来运行PhyML 3。对于250个序列的蛋白质模拟,这将PhyML的准确度提高到86.8%,这与FastTree的86.9%的准确度在统计上没有区别(,配对测试中,).我们还证实,FastTree的最小进化阶段产生的起始树比PhyML 3使用最大似然距离的BIONJ方法或RAxML实现简约(表1).

基于cat的分支长度和本地支持值

因为FastTree 2没有详尽地优化可能性,并且因为它报告的分支长度和局部支持值是使用CAT近似估计的,我们将其分支长度和局部支持值与基于长度和支持。具体来说,对于250个序列的蛋白质模拟,我们重新优化了分支长度,并使用PhyML 3和计算了FastTree拓扑的局部sh样支持值.(对于这两种工具,我们都使用了氨基酸进化的琼斯-泰勒-索顿(JTT)模型。)PhyML的内部分支长度与FastTree (= 0.90)。对于分支长度为1.0或更短的分支,平均差异仅为0.01,而对于分支长度在0.01和1.0之间的分支,平均差异为13%。对于正确分割的内部分支长度,FastTree与真实长度的一致性稍好(PhyML的中位数绝对差为0.062,FastTree为0.059)。因此,CAT近似给出了可接受的分支长度。

然而,如果精确的分支长度是必要的,那么CAT近似和标准近似是充分的。的在只有10个序列的比对中引入了近似,在更大的比对中,四个离散的速率类别可能不足以给出准确的可能性[15][22].对于16S核糖体rna的排列,分支长度可以短两倍分支长度(图S1).如在图S1,通过平均后验率进行校正可以减少这个问题,FastTree可以计算出快速但准确的近似的长度。

本地sh样支持度值在FastTree和PhyML (= 0.90)。对于FastTree或PhyML的局部支持值至少为0.9的分割,平均绝对差值仅为0.008,这并不比抽样误差大多少。(例如,在1000个自举和95%支持的情况下,由于抽样而产生的支持值的标准差为= 0.007)。FastTree在区分正确和不正确的分割方面不如PhyML有效,但差异很小:接受者操作曲线下的面积(AOC)为0.880,而不是0.887 (,测试[23]).

启发式的有效性

然后我们研究了FastTree 2的拓扑精度是如何受到其启发式的影响的。如表1例如,使用线性SPR的FastTree的最小进化阶段不如执行穷举SPR移动的最小进化方法FastME 2准确[8][12].FastME计算内部节点之间的距离与FastTree的最小进化阶段不同:FastME使用序列之间距离的平均值,而FastTree使用配置文件之间的距离,这是序列的平均值。然而,只有NNI移动的FastTree 1与只有NNI移动的FastME的结果非常相似[2].因此,我们将最小进化方法与SPRs在精度上的适度差异归因于FastTree的启发式。为了消除这种影响,我们使用FastME开始树运行FastTree。为了消除FastTree的ML启发式的影响,我们用穷举的ML NNIs运行它,并在每个NNI中对分支长度进行更穷举的优化(为每个四方优化分支长度4轮,而不是1-2轮)。结合FastME+SPR起始树和穷举NNIs, FastTree 2将5000个蛋白质序列的模拟比对精度从84.3%提高到85.0%。这种适度的影响说明了FastTree的所有启发式对准确性几乎没有影响,并且与添加ML SPRs相比,删除它们对拓扑结构的改善很少(例如,RAxML 7.2.1的准确率为88.4%)。

我们还对FastTree进行了超过78000个核苷酸序列的模拟测试。这些模拟来自16S核糖体RNA比对(见方法).这些模拟对齐的大尺寸使它们成为对FastTree启发式的严格测试。在这些模拟中,FastTree给出了比精确邻居连接或Clearcut更精确的拓扑[24],一种更快的启发式邻居连接变体(表1).(为了用精确的邻居连接分析如此大的对齐,我们使用了NINJA[25])。为了验证FastTree邻居连接阶段的启发式方法不会降低准确性,我们还以精确的邻居连接树作为起始树运行FastTree,然后进行最小进化NNIs、spr和ML NNIs。这提供了与常规FastTree相同的准确性,或者具有邻居连接阶段启发式最快设置的FastTree (-fastest)。这三个变体的拆分正确率为92.10%。

FastTree在不比较所有序列对的情况下也能达到精确的拓扑结构,这似乎令人惊讶。然而,最小进化NNIs和spr是“一致的”——它们能找到正确的树,即使距离包含一些误差,只要误差远远小于内部分支长度[26][27].在实践中,误差通常大于内部分支长度,但这仍然可能解释为什么NNIs和spr足以正确地找到大多数分割。

生物对齐拓扑的质量

为了确认FastTree为真正的排列找到了良好的拓扑结构,而不仅仅是在模拟中,我们在16S核糖体rna和COG的蛋白质家族上进行了测试。虽然这些家族相当大(分别有300,000或19,000个成员),但我们首先测试的只是500个序列的随机子集,以便我们可以运行PhyML 3.为了测量来自FastTree 2、PhyML 3和RAxML 7的拓扑结构的质量,我们使用a重新优化了分支长度模型(使用RAxML),并比较结果的可能性。从模拟中,FastTree发现了比PhyML 3更好的拓扑结构和NNI移动,但不如RAxML 7 (表2).

然后,我们在16S rrna和cog的更大比对上测试了FastTree和RAxML。对于数千个序列的对齐,RAxML 7.0.4有点慢,因此我们使用了RAxML 7.2.1,它引入了快速收敛选项以及其他优化。通过快速收敛,如果在一轮SPR移动中拆分的变化小于1%,RAxML将终止搜索。如表1对于5000种蛋白质,快速收敛的RAxML是相当准确的。

在较大的比对中,RAxML 7.2.1的可能性远高于FastTree,所有的可能性差异都具有统计学意义(所有, SH测试使用CONSEL[28]).然而,FastTree确实在RAxML拓扑中发现了大多数具有强大支持的分裂(表3).例如,FastTree发现96-98%的RAxML拆分的全局引导值为90%或更高。

缩略图
表3。比较RAxML和FastTree的对数可能性,以及FastTree与RAxML良好支持的分割的一致性,以获得真正的大型对齐。

https://doi.org/10.1371/journal.pone.0009490.t003

运行时间和内存要求

最后,我们比较了FastTree、RAxML和PhyML在真正对齐上的计算性能。如表4对于500个序列的比对,FastTree在使用相同的进化模型时比RAxML 7.0.4快100倍左右,相对于PhyML 3甚至更快。对于数千个序列的比对,FastTree仍然比RAxML 7.2.1快100-800倍,具有快速的SPRs收敛,而PhyML 3没有在合理的时间内完成。

对于目前最大的序列之一,包含237,882个16S核糖体rna, FastTree在台式电脑上花费了不到一天的时间和5.8 GB的内存。作为比较,假设RAxML仅花费2天时间处理15,011个序列,并且乐观地假设O()扩展,RAxML将需要大约半年的时间来完成完整的16S对齐。用传统的基于距离矩阵的最小进化方法分析这样的对齐也将是令人望而却步的——仅仅计算和存储这些序列的所有成对距离,而不计算拓扑,将需要大约一天半的时间和113 GB的存储空间。

所有的FastTree时间都包括本地sh样支持值的计算,而其他工具在没有支持值的情况下运行。本地支持值对FastTree的运行时间影响不大。例如,在7个包含2500个蛋白质序列的COG比对中,FastTree推断一棵树的平均时间是345秒,而计算sh样支持的平均时间是51秒。对于237,882个16S rrna的完全对齐,支架只需要一个小时。

的分支长度优化在RAxML 7.2.1中花费了大量的时间模型,即使CAT近似被用来搜索一个好的拓扑。(RAxML还可以在模型,但我们只在CAT模型下运行SPR移动,然后在CAT模型下优化分支长度,因为RAxML 7.2.1没有报告基于cat的分支长度。)如果不需要分支长度,例如在引导过程中,那么RAxML的速度可能比中所示的快2-3倍表4.例如,对于15,011个16S rrna,如果阶段被移除,那么RAxML 7.2.1需要30个小时而不是64个小时,这仍然比FastTree慢45倍。的RAxML阶段也预计将使所需的内存增加四倍。例如,对于15,011个16S rrna, FastTree需要0.56 GB内存,而RAxML需要需要2.6 GB。

随着时间的推移可能性的提高

为了更直接地比较FastTree和RAxML的搜索策略,我们比较了它们对4114个16S rrna核苷酸比对的可能性随时间的改善[11]对2500人的COG家族的7个蛋白质序列进行了分析。我们使用CAT近似和核苷酸取代的广义时间可逆(GTR)模型或氨基酸取代的JTT模型来运行这两种方法。我们用RAxML计算中间树和最终树的可能性,重新优化分支长度,以及图2显示了FastTree最小进化树和最终树、RAxML的初始简约树和连续几轮SPR移动的运行时间和对数似然,以及以FastTree最小进化树为起始树的RAxML的运行时间和对数似然。这些时间不包括FastTree的支持值或RAxML优化分支长度

缩略图
图2。随着时间的推移,真正的结盟的可能性。

每一行显示了使用不同工具达到给定可能性所需的时间。对于COG对齐,所有的时间和可能性都是七个对齐的平均值。对于FastTree,我们展示了最小演化拓扑和最终(近似ml)拓扑的时间和可能性的改进。对于RAxML,我们展示了最大简约开始拓扑、前两轮SPR移动和最终拓扑(注意中断的拓扑)轴)。对于具有FastTree(最小进化)开始树的RAxML,我们显示了开始拓扑和RAxML的前两轮SPR移动。

https://doi.org/10.1371/journal.pone.0009490.g002

给定相同的起始树,FastTree的ML阶段提高的可能性与一轮RAxML的SPR移动大致相同,并且在大约40%的时间内(图2).FastTree的ML阶段在寻找良好支持的分割(图S2).我们在其他大型16S对准中获得了类似的结果(表S1).虽然这个比较表明FastTree最初比RAxML快,但RAxML的第一轮SPR移动只占其运行时间的一小部分。FastTree和RAxML之间的速度差异主要是因为RAxML对更好的拓扑进行了更彻底的搜索,而RAxML则是分支长度。

起始树:最小进化vs最大简约

RAxML的简约阶段比FastTree的最小演化阶段慢4-17倍,通常比带有ML NNIs的FastTree慢。FastTree的速度优势随着更大的对齐而增长(数据未显示),这是预期的,因为FastTree应该按O()和RAxML的精简阶段使用随机的逐步相加,其规模为O(),以及有限的基于吝啬的SPR移动。还有更快的实现,比如RAxML 7.2.5(在我们进行上述实验之后发布)或TNT[29],但这些仍然按O().对于15011个16s的rna, RAxML 7.2.5的精简和FastTree的最小进化阶段大约需要相同的时间(数据未显示)。

根据可能性的测量,FastTree的最小进化开始树在COG比对中比RAxML的简约开始树好得多,但在大型16S rRNA比对中要差得多(图2而且表S1).可能性的差异反映了标准,而不仅仅是搜索策略的差异:对于COG对齐,RAxML简约起始树比FastTree的最小进化树更简约(平均简约得分分别为281,237和283,125)。相反,对于4114个序列的16S比对,FastTree的最小进化树比简约树短(长度分别为43.0和44.6)。对于这种对齐,最小进化树的对数可能性比简约树差2705,但最小进化树在最终的RAxML树中发现了更多强支持的分割:最小进化树通过全局引导发现了851个分割中的826个90%,而吝啬发现了814个。因此,我们不确定这种可能性的差异在生物学上是否有意义。

讨论

我们已经证明,FastTree 2可以在合理的时间内计算出精确的拓扑结构,对多达数十万个序列进行比对。FastTree是开源软件,可以在http://microbesonline.org/fasttree.C源代码有广泛的文档记录,欢迎您的贡献。每个微生物基因家族的FastTree树,包括有数万个成员的家族,如ABC转运体,可在MicrobesOnline (http://microbesonline.org/),以及用于检查这些树的“树浏览器”。这些树将在MicrobesOnline的下一个版本中从FastTree 1更新到FastTree 2。

由于DNA测序技术正在迅速发展,我们预计很快就能对数百万个序列进行比对。对于这些巨大的对准,最需要计算的步骤将是最初的邻居连接阶段。在FastTree 2.0中,邻居连接需要O()时间和O()空间,而其他阶段最多占用O()时间和O()空间。例如,对于237,882个16S序列,FastTree 2.0的邻居连接阶段已经花费了21.8小时中的10.8小时。在FastTree 2.1中,我们改进了时间和内存的比例,从O()至O(),而不会影响模拟的准确性(数据未显示)。FastTree 2.1还支持并行执行邻居连接阶段的关键步骤。为了进一步提高可伸缩性,有可能使用分治方法在O()时间,如在PartTree[30].在我们的模拟中,PartTree起始树不允许FastTree达到与FastTree的邻居连接起始树相同的精度(数据未显示),但是分治方法可能仍然足以获得部分解析的初始树。

如此庞大的家族也为多序列比对提出了挑战。我们使用剖面比对来避免在大家族上进行多序列比对的挑战。这对于16S rna很有效,因为Infernal利用了高度保守的二级结构[1],但我们不确定它对不同的蛋白质家族给出了准确的结果。相比之下,传统的渐进式多序列对齐方法是不可扩展的,因为它们的输出增长为O():有O()独立的插入,并且每次插入都需要在对齐中添加一个新列,因此O()存储。然而,快速统计对齐使用O()表示,包括内部和作为输出格式[31].将这种表示与快速引导树结构相结合,应该可以用数百万个序列构建渐进的多个序列对齐。

最后,如何评估如此大的树木的质量或可靠性尚不清楚。不同的方法给出了非常不同的拓扑结构和很大的可能性差异,但很少有差异得到了bootstrap的良好支持。事实上,一个可能性相对较差的拓扑结构可能与最佳树有相对较好的一致性。这可能表明,可能性较高的树包含许多改进,但很少有个别改进具有统计学意义。如果系统发育信号有限,这是预期的。或者,bootstrap可能过于保守。本地支持值确实显示出更多的显著差异(表3),但局部支持值可能会向上偏置,因为它们没有考虑所有的替代拓扑。这些问题需要进一步研究。

材料与方法

最小进化“线性”子树修剪-再嫁接

为了减少从O()至O(), FastTree只做两轮“线性SPRs”。对于每个节点,FastTree会穷尽搜索长度为2的移动。它将每一步移动的长度都扩展到10沿着沿途的每个点上的最佳选择。

正如Richard Desper和Olivier Gascuel提出的,FastTree将每个潜在的SPR移动视为NNIs序列。SPR移动的树长度变化只是NNIs变化的总和,就像

NNI的树长度变化,其中A, B, C, D可能是子树而不是序列,由.在FastME中,引入了平衡最小进化[12]是a和c的成员之间距离的拓扑加权平均值。相反,在FastTree中,子树A和C的轮廓线与轮廓线之间的对数修正距离是多少子树的子树是由子树的子树派生而来的.(虽然FastTree 1使用加权连接,就像在BIONJ中一样,FastTree 2使用非加权连接,因为它们更快,而且对精度的轻微影响被ML NNIs消除了。)对于核苷酸序列,对数修正是Jukes-Cantor修正,在那里是各剖面位置的平均不相似度。对于氨基酸序列,FastTree使用类似于记分器的经验对数校正[32],在那里是基于从BLOSUM45相似矩阵衍生出来的氨基酸不相似矩阵。

在FastME中,上面关于树长度变化的公式是完全正确的,因为树中其他分支长度的变化可以表示为相互抵消的距离的组合[26].然而,在FastTree中,树长度变化的公式是一个近似值,因为对数修正的距离不会以这种方式抵消。然而,带有NNIs的FastTree和带有NNIs的FastME给出了非常相似的结果[2],并且计算总树长度的确切变化并不能提高FastTree的SPRs的准确性(数据未显示)。

最大似然阶段

最大似然相的关键数据结构是树拓扑结构、分支长度和每个内部节点的后验分布。(FastTree存储的树在根处有一个三叉分支,但根的位置在生物学上没有意义,也不影响可能性[5])。给定分支长度和下面的序列,内部节点的后验分布描述了相应祖先节点的状态。例如,对于核苷酸数据,它存储给定位点是a、C、G或t的概率。FastTree存储后验分布内部节点(不是根节点),它们需要O()空格,其中对齐的长度和是字母表中的字符数。

关键的原语操作是(1)计算两个后验分布的联合似然,给定它们之间的长度,(2)计算父节点的后验分布,给定它的两个子节点的后验分布和它们的两个分支长度。这些足以计算出树的可能性[5]:如树(A,B,(C,D))的似然为其中AB和CD为后验分布。

在ML阶段的开始,我们有一个最小演化拓扑和分支长度。最大似然阶段的步骤如下:

  • 计算每个节点的近似后验分布,使用其子节点的加权平均值。虽然最初的后验分布是近似的,但所有未来对拓扑或分支长度的更改都会将后验分布更新为它们的确切值。
  • 优化一个回合的所有分支长度,使用不带参数的简化模型(不带CAT,如果需要GTR,则使用Jukes-Cantor而不是GTR)。
  • 使用简化模型执行一轮ML NNIs。
  • 如果正在使用GTR模型,则优化核苷酸跃迁速率参数,从Jukes Cantor模型切换到GTR模型并重新计算后测分布,并使用新模型优化一轮的所有分支长度。
  • 如果正在使用CAT模型,请估计每个站点的费率类别,重新计算后验分布,并使用新模型优化一轮的所有分支长度。
  • 执行额外的ML NNIs轮,包括子树跳过和星形拓扑测试。
  • 执行最后一轮ML NNIs,不跳过子树或星型拓扑测试。
  • 优化一回合的所有分支长度。
  • 计算类sh的本地支持值。

一轮ML NNIs。

在每一轮NNIs中,FastTree在访问父节点之前访问每个节点(深度优先后序遍历)。在每个节点上,它比较树的可能性,,其中A和B是它的子结点,C是它的兄弟结点,D是树的其余结点。在这个过程中,FastTree使用A、B和C的后验分布(或序列)和“向上后验”D,它表示树的其余部分。更精确地说,D是给定所有节点的父节点N的后验分布N的孩子(见图3).这些向上可以被认为是一种暂时在当前位置重新根树的方法。特别地,树的似然值可以从后面的A、B、C和D计算出来。

缩略图
图3。穿越一棵向上的树。

FastTree通过分析子树A、B和C的后验分布以及“向上-后验”D来优化节点N附近的树。

https://doi.org/10.1371/journal.pone.0009490.g003

一个节点的后验分布可以通过其父节点的后验分布和其兄弟节点的后验分布来计算。FastTree只存储这些从树中当前位置到根的路径的向上路径,所以它们取O()空间,其中是树的最大深度。因为FastTree总是在父节点访问子节点之前访问子节点,所以它使用的后验分布和上后验分布是最新的,即使在拓扑结构发生变化时也是如此。

当它访问每个节点时,对于每个节点周围的三个备选拓扑,FastTree会优化分支长度以使可能性最大化。对于拓扑,从当前树中设置5个初始分支长度。对于其他拓扑,保持到A、B、C和D的分支长度,就像内部分支长度一样。假设有一个四重奏(比如说), FastTree首先优化AB和CD之间的分支长度,然后优化通往A、B、C和d的分支长度。FastTree优化每个分支长度的精度为0.0001或0.1%,以较大者为准。这五个优化定义了四重奏的一轮优化。在一轮优化中,FastTree重用了一些内部后验分布:它需要AB和CD的后验分布,以便优化AB和CD之间的分支长度,然后它需要BCD的后验分布,ACD, AB的新后验分布,给定a和B的新分支长度,最后是ABD和ABC。

默认情况下,FastTree在一个回合中优化所有三个四重奏拓扑中的分支长度。任何明显比当前拓扑差(5个对数似然单位)的拓扑都在第一轮之后被放弃。如果存在多个拓扑,则对剩余的拓扑进行下一轮优化。几轮优化完成后,FastTree会在必要时更新拓扑。在任何一种情况下,它都会将分支长度更新为重新优化的值,并重新计算节点的后验分布。

log概率5的差异看起来可能是一个很小的差异,因此启发式可能会错过对拓扑结构的一个很好的更改。然而,在第一轮之后对分支长度进行优化通常会导致对数似然的小幅改善。例如,如果我们使用FastTree和GTR+CAT模型分析40个随机选择的16S rrna,并将分支长度优化轮增加到4轮(-mlacc 4),那么任何NNI在第二轮分支长度优化中的平均改进仅为1.1个对数似然单位,在第3轮和第4轮组合中仅为0.04。为了正确地看待这些数字,小于2的对数似然差异在统计上并不显著((似然比检验),常见似然变化大得多的NNIs。对于含有5000个蛋白质序列的模拟比对,始终优化两轮可以使准确度提高微不足道(0.03%),并使运行时间增加23%。

优化模型参数。

在第一轮NNIs之后,FastTree优化模型中的任何参数。首先,如果使用GTR模型,则有六种相对速率需要优化,每种核苷酸转换一个。(过渡矩阵的平稳分布被设置为四个核苷酸的经验频率。)FastTree通过依次数值优化模型中的六个参数来优化树的可能性(具有固定的分支长度和拓扑结构)。随着模型中的每一个变化,它重新计算所有的后验分布。然后它会第二次优化这六个参数。这并没有完全优化模型参数,但它给出了可接受的结果(表2).

其次,除非设置了-nocat选项,否则FastTree会估计每个站点的进化速度。给定所需种类的相对速率FastTree选择值之间的对数间隔而且.默认情况下,,相对比率在0.05 ~ 20之间。对于这些相对率,FastTree重新计算所有后验分布,并重新计算树在每个站点的对数似然值。FastTree然后使用贝叶斯方法来选择在每个站点使用的速率:FastTree最大化,在那里是分布先验。为了避免过拟合,我们使先验在对齐中比实际率变化更有峰值:先验的形状参数为3,尺度参数为1/3,平均值为1。在选择了费率类别后,FastTree会调整费率,使所有站点的平均费率为1.0。

我们确认,贝叶斯方法设置的速率类别防止过拟合在小对齐。例如,在只有10个序列的模拟蛋白质比对中(从[2]),加入CAT模型后,FastTree的准确率由76.2%提高到78.0%。(为了进行比较,PhyML没有或SPRs的准确率为74.4%[2])。相反,在核苷酸模拟中,24个序列(不现实地)不包含任何跨位点的速率变化(快速进化的排列)[12])时,CAT模型的准确率仅略有下降,从93.6%下降到93.4%。(为了进行比较,PhyML没有或SPRs的准确率为93.6%[2])。

完成ML NNIs。

在后来的NNIs中,FastTree使用了更精确的模型,并使用了两个额外的启发式方法“子树跳过”和“星型拓扑测试”,下面将对此进行描述。如在结果,这些启发式对准确性的影响很小。

如果没有NNI导致任何四项的可能性提高超过0.1,那么FastTree认为NNI已经收敛。FastTree重复NNIs轮直到收敛,直到达到轮次,需要O()时间。这是ML阶段最慢的部分。轮数的限制确保了可预测的运行时间,但FastTree通常会在达到限制之前收敛,即使是对于237,882个16S rRNA序列这样的大型对齐。我们选择了限制,这样错位的子树可以在(大致平衡的)树中移动,因子2是一个任意的安全因子。

收敛后,FastTree执行最后一轮ML NNIs,与第一轮一样,跳过子树并关闭星型拓扑测试。我们认为这是启发式的安全阀。最后,FastTree进行最后一轮优化分支长度,并计算类sh的本地支持。

子树跳过。

子树跳过背后的直觉是,如果子树在最近几轮NNIs中没有改变,那么进一步优化子树的尝试将是徒劳的。具体来说,在ML NNIs期间,FastTree不会遍历在前两轮中都没有看到任何显着的似然改善(0.1 log似然单位)的子树。在跳过子树之前,FastTree还会检查与父节点相邻的节点是否在上一轮中受到显著改进的NNI的影响。“子树跳过”启发式通常提供3倍的加速,使其成为FastTree的ML启发式中最重要的。子树跳过可能对SPR移动也有用。

星型拓扑测试。

如果当前拓扑(A,B,(C,D))比星型拓扑(A,B,C,D)好得多,那么NNI不太可能给出改进。具体来说,如果当前拓扑结构(优化内部分支长度后)明显(5个对数似能单位)比星型拓扑结构更有可能,那么FastTree不会优化其他分支长度或考虑两种备选拓扑结构。但是,FastTree只在节点在最后一轮NNIs中没有改变时才使用这种启发式方法。为了近似星形拓扑的似然,FastTree使用内部分支长度最小为0.0001的似然。

分支长度。

为了在ML阶段的开始和结束以及优化模型参数之后优化树中的所有分支长度,FastTree再次使用后序遍历。在每个节点上,它考虑节点的子节点和父节点上的三个节点星形拓扑,对两个子节点使用后验分布,对自身使用上-后验分布。(在根节点上,它使用了所有三个子节点。)它在数值上优化了这三个分支长度连续两轮。

类似sh的本地支持。

对于每个节点,本地支持来自当前拓扑和两个备用(NNI)拓扑的每个站点可能性。对于当前拓扑,FastTree使用当前(已经优化过的)分支长度。对于替代拓扑,FastTree为四重奏优化分支长度,就像在NNIs期间一样,最多可以优化两轮。给定这三种拓扑的每个站点可能性,FastTree使用具有1,000个自举重复的SH检验来估计给定分割的置信度[19].如果附近存在分辨率较低的节点,则应该谨慎地解释支持值,因为可能没有考虑高可能性的替代拓扑。

可能性计算的低级优化。

RAxML存储似然向量(即子树和内部节点上给定字符的联合似然),FastTree存储后验分布,这些分布经过标准化处理,使每个站点的值之和为1。这可能会提高大型对准的数值稳定性。为了减少内存使用,FastTree将这些向量存储在单精度浮点数中。树或特定地点的对数可能性以双重精度存储。

与RAxML类似,FastTree以旋转的形式存储后验分布,并乘以转换矩阵的特征矩阵。(对于Jukes Cantor模型,这是不必要的。)这减少了从O()至O(),而将计算后验分布的代价留在O()(但常数系数较高)。

在计算一对后验分布的联合似然时,FastTree通过操作似然而不是对数似然来避免在每个站点执行对数。为了防止数值下溢,FastTree在必要时通过常数重新调整可能性。每当执行此操作时,它都会更新一个单独的(基于对数可能性的)计数器。类似地,当在每个站点上计算树的可能性时,例如在优化费率类别时,FastTree在访问每个节点后,如有必要,会重新调整每个站点的可能性。

FastTree使用SSE2指令(Intel和AMD最近的cpu的一个特殊功能),用一条指令操作4个单精度浮点值。这将加快蛋白质校准计算速度高达50%(数据未显示)。

数值优化。

为了找到优化可能性的参数,FastTree使用Brent的方法,这是一种数值方法,迭代地将它正在搜索的区间减半(http://en.wikipedia.org/wiki/Brents_method).因为Brent的方法只在一维中操作,FastTree依次优化不同的参数,然后重复优化(例如,它优化第一个分支长度,然后是第二个,然后是第三个,然后重复)。

FastTree根据最初的猜测估计搜索的初始间隔(例如,分支之前的长度)和交替值而且.如果小于最小值,它使用最小值,,代替。如果初始猜测没有包含最小值(也就是说,中间值并不比两个端点好),那么FastTree将扩展搜索间隔,直到包含最小值为止。然而,小的间隔通常是足够的。FastTree还会在参数变化很小或比例很小的情况下终止优化。总之,这些修改消除了大约三分之一的可能性评估。

生物和模拟对准

模拟的蛋白质排列和真实的COG排列已经在前面描述过[2].237,882个不同序列的16S比对来自GreenGenes[33]http://greengenes.lbl.gov).与15,011个不同“家族”的16S比对是这些序列的非冗余子集(相同的)。与500个序列的16S对齐也是非冗余随机子集(相同的)。其他大型16S对准是从[11]

对于具有78132个不同序列的16S类模拟,我们使用了从完整16S序列集的非冗余对齐子集推断出的最大似然树(使用Jukes-Cantor模型(没有CAT)的FastTree(1.9)的早期版本。为了确保模拟树是可解析的,这有助于方法的比较(但会提高所有方法的准确性),小于0.001的分支长度被替换为0.001的值,这大致相当于在内部分支上进行一次替换,因为16S对齐有1287个位置。根据变异系数为0.7的gamma分布,从16个速率类别中随机选择每个站点的进化速率。给定树和速率,用Rose模拟序列[34]在HKY模型下,不存在转移偏差。为了允许Rose处理小于1%的分支长度,我们设置“MeanSubstitution = 0.00134”并将分支长度乘以1000。

软件使用

我们使用FastTree 2.0.0。我们使用了2009年7月6日发布的PhyML 3.0源代码,并从1修改了BL_MIN。E-10到1。e-8来克服一些模拟蛋白质对齐的数值问题,如Stepháne Guindon所建议的。FastME 2.06由Olivier Gascuel提供。RAxML 7.0.4和7.2.1从作者的网站获得。RAxML 7.2.1使用SSE指令编译。NINJA由Travis Wheeler提供,可在http://nimbletwist.com/software/ninja/.BIONJ由http://www.lirmm.fr/~w3ifa/MAAS/BIONJ/BIONJ.c.BIONJ使用与phylip's protdist (http://evolution.genetics.washington.edu/phylip.htm)和JTT模型(无gamma)。使用FastTree和-makematrix选项获得对数校正距离。

支持信息

图S1。

200个16S rRNA序列比对的分支长度随所使用的Γ近似值而系统地变化。CAT长度来自FastTree,所有Γ分支长度来自PhyML,具有FastTree的拓扑结构和优化的形状参数。顶部面板显示,来自各种模型的分支长度彼此之间有一个大致的线性关系,但它们有不同的尺度。底部面板显示了树的总长度如何随着类别数量的变化(注意log χ轴)。“Use Median”长度来自使用-use_median运行PhyML,它使用每个区域的中值而不是平均值来近似gamma分布。“校正”长度是“使用中位数”长度乘以平均后验率,可以通过使用-print_site_lnl运行PhyML来获得(感谢Stepháne Guindon指出这一点)。修正后的长度比其他速率更快地收敛到正确的值。“CAT/Gamma”树长度,来自FastTree 2.1的-gamma,也是相当准确的。通过这个选项,FastTree 2.1优化了Γ20.使用基于CAT模型优化的FastTree的20个相对速率和分支长度的站点可能性。

https://doi.org/10.1371/journal.pone.0009490.s001

(0.13 MB ps)

图S2。

总的分割或强支持的分割与RAxML的Final树不一致,相对于时间。16S树有4111个分支,而COG树各有2497个分支。COG树的所有值都是7个COG的平均值。

https://doi.org/10.1371/journal.pone.0009490.s002

(0.02 MB ps)

表S1。

16S rRNA大排列的次数和可能性

https://doi.org/10.1371/journal.pone.0009490.s003

(0.02 MB pdf)

致谢

我们感谢Alexandros Stamatakis提出了FastTree和RAxML之间的时间和可能性比较,并对手稿进行了评论。

作者的贡献

构思并设计实验:MNP PD。进行实验:MNP。分析数据:MNP。撰写论文:MNP PD APA。

参考文献

  1. 1.Nawrocki EP, Kolbe DL, Eddy SR (2009) Infernal 1.0: RNA排列的推断。生物信息学15:1335-7。
  2. 2.Price MN, Dehal PS, Arkin AP (2009) FastTree:用轮廓代替距离矩阵计算大型最小进化树。Mol Biol vol 26: 1641-50。
  3. 3.齐藤,内敏(1987)邻接法:一种重建系统发生树的新方法。Mol Biol Evol 4: 406-425。
  4. 4.研究者JA, Keppler KJ(1988)关于Saitou和Nei的邻居连接算法的注释。Mol Biol Evol 5: 729-31。
  5. 5.Felsenstein J (1981) dna序列进化树:最大似然方法。中华药理学杂志17:368-376。
  6. 6.用最大似然重建系统发生树是困难的。IEEE/ACM计算机计算生物学生物信息学3:92-94。
  7. 7.Guindon S, Gascuel O(2003)一种通过最大似然估计大系统发育的简单、快速和准确的算法。系统生物学52:696-704。
  8. 8.Hordijk W, Gascuel O(2005)提高基于最大似然的系统发生树搜索算法中SPR移动的效率。生物信息学21:4338-4347。
  9. 9.Stamatakis A (2006) RAxML-VI-HPC:基于最大似然的系统发育分析,涉及数千个类群和混合模型。生物信息学22:2688-2690。
  10. 10.Felsenstein J(1985)系统发育的置信限制:一种使用自举法的方法。进化39:783-791。
  11. 11.Stamatakis A, Hoover P, Rougemont J (2008) RAxML web服务器的快速引导算法。系统生物学57:758-771。
  12. 12.Desper R, Gascuel O(2002)基于最小进化原理的快速准确的系统发育重建算法。计算生物学杂志9:687-705。
  13. 13.Nei M, Kumar S, Takahashi K(1998)在系统发育分析中,当使用的核苷酸或氨基酸数量较少时,优化原则往往会给出不正确的拓扑结构。美国自然科学杂志95:1290 - 7。
  14. 14.速率异质性的系统发育模型:一个高性能计算的角度。
  15. 15.杨智(1994)基于可变比例DNA序列的最大似然系统发育估计:近似方法。中华药理学杂志39:306-314。
  16. 16.金登,戴修,杜法亚德,加斯库尔O(2009)利用PhyML估计最大似然系统发育。方法Mol Biol 537: 113-37。
  17. 17.Guindon S, Dufayard J, Lefort V, MAnisimova, Hordijk W等(2010)估计最大似然系统发生的新算法和方法:评估PhyML 3.0的性能。系统医学杂志。:在印刷中。
  18. 18.Gascuel O (1997) BIONJ:基于序列数据简单模型的NJ算法的改进版本。Mol Biol Evol 14: 685-695。
  19. 19.Shimodaira H, Hasegawa M(1999)对数可能性的多重比较及其在系统发育推断中的应用。Mol生物进化16:1114-1116。
  20. 20.Tatusov RL, Natale DA, Garkavtsev IV, Tatusova TA, Shankavaram UT,等(2001)COG数据库:全基因组蛋白质系统发育分类的新进展。核酸Res 29: 22-8。
  21. 21.Anisimova M, Gascuel O(2006)分支的近似似然比检验:一种快速、准确和强大的替代方法。系统生物学55:539-52。
  22. 22.马可夫链与分子进化的协变过程。计算生物学报11:727-733。
  23. 23.DeLong ER, clark - pearson DL(1998)比较两条或多条相关的受试者工作特征曲线下的面积:一种非参数方法。生物识别技术44:837-45。
  24. 24.刘建民,李建民,李建民,等(2006)一种基于距离的快速系统发生树构建方法。中华药理学杂志62:785-92。
  25. 25.Wheeler TJ(2009)使用NINJA进行大规模邻居连接。
  26. 26.Desper R, Gascuel O(2004)系统发育推断平衡最小进化方法的理论基础及其与加权最小二乘树拟合的关系。Mol Biol Evol 21: 587-598。
  27. 27.bordwich M, Gascuel O, Huber KT, Moulton V(2009)基于系统发育推断平衡最小进化原理的拓扑运动一致性。IEEE/ACM计算机计算机生物信息学学报6:110-7。
  28. 28.张志刚,张志刚(2001)系统发育树选择的可信度评价。生物信息学17:1246-1247。
  29. 29.Goloboff PA, Farris JS, Nixon KC (2008) TNT,一个免费的系统发育分析程序。分支学24:774-786。
  30. 30.Katoh K, Toh H (2007) PartTree:一种从大量未对齐序列构建近似树的算法。生物信息学23:372-374。
  31. 31.Bradley RK, Roberts A, Smoot M, Juvekar S, Do ea J(2009)快速统计对齐。PLoS计算生物学5:
  32. 32.Sonnhammer ELL, Hollich V (2005) Scoredist:一个简单而健壮的蛋白质序列距离估计器。BMC生物信息学6:108。
  33. 33.DeSantis TZ, Hugenholtz P, Larsen N, Rojas M, Brodie EL等(2006)Greengenes,嵌合体检查的16S rRNA基因数据库和与ARB兼容的工作平台。应用环境微生物72:5069-5072。
  34. 34.王志刚,王志刚,王志刚(1998)序列族的生成。生物信息学14:157-163。