[文献分享]Molecular phylogenetics: p

首段放文章链接：https://www.nature.com/articles/nrg3186

最近在帮一个师姐做直系同源基因的物种进化树，补充一些进化方面的基础知识是很有必要的，毕竟没有生物学背景的分析特别危险（https://mp.weixin.qq.com/s/5HXUH4lWVc70d1CKJOQ73g）。

Abstract

这篇文章主要总结了进化分析的一些主要方法，包括parsimony,distance,likelihood,和Bayesian methods。并且讨论了它们各自的优点和缺点，并对于它们的使用提供了指导。

Noun explanation

systematics(系统学):推断物种间的进化关系，并且使用这些信息来进行分类。

Taxonomy(分类学)：描述、分类和命名物种。

Coalescent(合并)： 追溯从现代人口中随机抽取的一个序列样本的系谱关系时，将祖先谱系连接起来的过程。

Gene trees(基因树)：在一个基因位点或基因组区域序列的系统发育或系谱树。

Statistical phylogeography（统计系统地理学）： 对密切相关物种的种群数据进行统计分析，以推断种群参数和过程，如种群大小、人口统计学、迁移模式和速率。

Species tree（物种树）：基于单个位点基因树的一系列物种进化树。

Systematic errors（系统误差）： 由不正确的模型假设引起的错误。当数据量增加时，误差会加剧 。

Random sampling errors（随机样本错误）： 由于数据有限，参数估计中的误差或不确定性。

Cluster algorithm（聚类算法）： 一种将一组个体分配给组(或集群)的算法，使同一集群中的对象之间比来自不同集群的对象更相似。层次聚类分析可以是聚集性的(从单个元素开始，依次将它们加入集群)，也可以是分裂性的(从所有对象开始，依次将它们划分为分区) 。

Markov chain（马尔科夫链）： 一个随机的状态序列(或链)，其性质是，给定当前状态，下一个状态的概率不依赖于过去的状态。

Transitions（转换）：嘧啶之间(T↔C)或嘌呤之间(A↔G)的替换。

Transversions（颠换）：嘧啶和嘌呤之间替换(T or C↔A or G).

Unrooted trees（无根树）: 进化树中根的位置是不确定的。

Long-branch attraction: 一种不正确的树，其长枝干基于parsimony或基于过于简单的模型不正确地组合在一起。

likelihood ratio test: 使用likelihood来比较两个嵌套假设的一般的假设检验方法,通常使用χ2分布来评估。

Molecular clock: 一种假说，假定进化率是随着时间或跨系谱而保持不变。

Prior distribution(先验分布)： 在分析数据之前分配给参数的分布。

Posterior distribution(后验分布)： 参数(或模型)的分布取决于数据。它结合了先验和数据(likelihood)中的信息 。

Markov chain Monte Carlo algorithms(MCMC 算法)： 蒙特卡罗模拟是使用随机数对生物过程进行的计算机模拟。MCMC算法是一种蒙特卡罗模拟算法，它从目标分布(通常是贝叶斯后验分布)中生成样本。

直系同源，旁系同源的区别(From Wikipedia)

在遗传学中，同源这一概念主要是指序列同源，表明两个或多个蛋白质或DNA序列具有相同的祖先。同源的序列也很可能有相似的功能。 蛋白质和DNA的同源性常常通过它们序列的相似性来判定。

同源序列可分为两种：直系同源（orthology）和旁系同源（paralogy）。直系同源的序列因物种形成（speciation）而被区分开（separated）：若一个基因原先存在于某个物种，而该物种分化为了两个物种，那么新物种中的基因是直系同源的；旁系同源的序列因基因复制（gene duplication）而被区分开（separated）：若生物体中的某个基因被复制了，那么两个副本序列就是旁系同源的。直系同源的一对序列称为直系同源体（orthologs），旁系同源的一对序列称为旁系同源体（paralogs）。若两个旁系同源基因（paralogs）的基因拷贝（gene duplication）是已发生在物种形成前，则这两个旁系同源（paralogs）可互称为“前重复同源基因”（outparalogs），可简称为“前重复”（outparalogs）；若于物种形成后某个直系同源基因（an ortholog）才发生基因复制（gene duplication），则互称为“后重复”（inparalogs）。

直系同源体通常有相同或相似的功能，但对旁系同源体则不一定：由于缺乏原始的自然选择的力量，繁殖出的基因副本可以自由的变异并获得新的功能。

Introduction

分子系统进化学对于基因组比较来说，成为了一个不可或缺的工具。在本文中被用来对宏基因组序列进行分类；在新测序的基因组中鉴定基因、调节原件和非编码RNA；对现代和祖先的个体基因组进行解读；构建祖先基因组。

Phylogenetic tree reconstruction: basic concepts

每一个分支代表一个遗传谱系随时间的延续，每一个节点代表一个新谱系的诞生。 如果树表示一组物种之间的关系，则节点表示物种形成事件。 在从一个群体中取样的序列的基因树中，节点代表个体从祖先中的出生事件，而在一个旁系基因家族（paralogous gene families）的树中，节点可能代表基因复制事件。

进化树的构建要么是基于距离的，要么是基于字符的。 在基于距离的构建中，计算出每对序列之间的距离，并将得到的距离矩阵用于树的重建。例如，neighbour joining 应用聚类算法到距离矩阵来得到一个完全解决的种系发生事件。在基于字符的构建方法中，包maximum parsimony, maximum likelihoodand Bayesian inference methods。这些方法同时比对所有的序列， 每次只考虑一个字符(相应对齐的位置)来计算每棵树的分数。 “树分数”是maximum parsimony的最小变化数，maximum likelihood的对数似然值和贝叶斯推理的后验概率。理论上来说，通过比对所有可能的数能得到一个分数最高的树。实际上，除了小数据集外，大型数据是没法这样进行的，因此就采用了启发式的树搜索算法。一个启发式的树搜索并不能保证在准则下得到一个最好的树，但是它对于分析大数据集是很有帮助的。 为了描述数据，distance matrix、maximum likelihood和贝叶斯推理都使用了替代模型，因此都是基于模型的，而maximum parsimony没有一个明确的模型，其假设是隐含的。

001.png

上图解读：

系统发生是系谱史演化的一种模式，其中分支的长度是未知的参数。a图有两个node，代表了两次物种形成事件，分别发生在t0期和t1期。 分支长度(b0、b1、b2和b3)通常表示为每个位点的期望替换数量的单位，并沿着分支来衡量演化量。

如果取代率随时间推移或在谱系间的变化是恒定的，我们就说分子钟保持不变。这样的话，树因此就有根， 这意味着从树的顶端到根的距离都是相等的，也就是图上的（b0 + b1 = b0 + b2 = b3). 一个有s个物种的有根树可以用s - 1个祖先节点的年龄来表示，从而涉及到s - 1个分支。 通过假设分子钟来推断有根树的过程叫做分子钟有根树。 对于亲缘关系较远的物种，不应假设分子钟假说。大多数进化树并不使用分子钟假说。## 如果树上的每个分支都被允许有一个独立的进化率，那么常用的模型和方法就无法识别根的位置，因此只能推断出没有根的树。对于s个物种的无根树有2s-3个分支数。 一种常用的“使树生根”的策略是在分析中包含外群体物种，这些外群体物种的亲缘关系比相关物种的亲缘关系要远得多， 尽管所有物种的推断树是无根的，但根被认为位于通向外群的分支上，因此内群物种的树是有根的。这种策略被称为外群生根。

Distance matrix method

D.1 Distance calculation

以核苷酸取代的马尔科夫链模型为基础，计算成对的序列距离。JC69模型假定任意两核苷酸的替换是相同的比率，而K80模型假定转换(transitions)和颠换（transversions）是不同的比率。这两个模型都预测了四种核苷酸等同的频率。在HKY85模型和general time reversible(GTR)模型中，放宽了碱基频率相等的假设。 由于局部突变率和选择性约束的变化，DNA或蛋白质序列的不同位点往往以不同的速率进化。 在距离计算中,这样的速率变化是通过假设gamma(Γ)位点的分布率,导致模型如JC69 +Γ，HKY85 +Γ或GTR +Γ 。

002.png

D.2 Distance matrix methods

计算完距离后，在基于距离的算法中不再使用序列比对。这里我们提到了三种方法：least squares(最小平方), minimum evolution and neighbour joining。

least squares最小化了距离矩阵中树上的计算距离(dij)和预期的距离(dˆij)差异(也就是说,在树上分支长度的总和连接了这两个物种i和j)，公式为：

003.png

minimum evolution方法使用树的长度(即分支长度的总和)而不是Q来进行树的选择，即使分支长度仍然可以使用least squares来估计。在minimum evolution法则中，较短的树似乎比较长的树更可靠。

neighbour joining 是使用最广泛的计算距离的方法。这是一个cluster算法，从一个星形树开始，依次选择一对分类群（基于距离的分类群）连接在一起，直到得到一个完全解析的树。 选择加入的taxa是为了最小化对树长度的估计 。在MEGA中实施了有效的 neighbour joining。

004.png

D.3 Strengths and weaknesses of distance methods

基于距离算法的一大优势（特别是neighbour joining）是它们的计算效率。 聚类算法是快速的，因为它不需要和maximum parsimony、maximum likelihood那样，比较最优条件下的树的数量。neighbour joining对于分析具有低水平序列差异的大型数据集是有用的。 请注意，使用现实的替代模型来计算两两距离可能很重要 。 distance methods对于差异大的序列表现很差，因为大的距离包含了很大的样本误差，而且大多数距离方法(例如neighbour joining)不能解释大距离估计的高水平的变化。distance methods也对序列比对中的gaps十分敏感。

Maximum parsimony

M.1 Parsimony tree score

在进化树中，Maximum parsimony通过将字符状态分配给系统发生树的内部节点来最小化系统发生树的变化数量。 字符(或位点)长度是该位点所需的最小更改数，而tree score是所有位点上字符长度的总和。maximum parsimony树是最小化树分数的树。

通过parsimony进行的树的比较，其中一些位点并不是有用的。 例如，在所有物种中都出现相同核苷酸的恒定位点，在任何树上的字符长度都为零。 单例位点（其中只有一个物种有不同的核苷酸，而其他所有物种都是相同的），也可以忽略，因为字符长度总是1 。 parsimony信息量位点是指至少观察到两个不同的特征，每个特征至少两次。对于四个物种来说，只有三个位点类型提供了有用信息：xxyy, xyxy, xyyx,其中x和y是任一两个不同的核苷酸。对于四物种来说，有三种可能的无根树， 哪一种类型是Maximum parsimony，取决于三个位点模式中的哪一个在比对过程中最常发生。

Fitch和Hartigan开发了一种算法，用于找出二叉树的最小变化数(以及重建原始状态以达到最小) 。PAUP,MEGA,TNT是通常使用的parsimony 程序。

Maximum parsimony的使用仍然很普遍:不是因为它被认为是没有假设的，而是因为它经常产生合理的结果，并且计算效率很高。

M.2 Strengths and weaknesses of parsimony

parsimony 的一大优点是它的简约性，它容易去描述和理解，并且它服从严格的数学分析。

parsimony的一个主要缺点是缺乏明确的假设，这使得在树的重建中几乎不可能包含任何序列进化过程的知识 由于在同一地点未能纠正多个取代物的错误，由此导致了一个被称为long-branch attraction的问题。如果正确的树由一个短的内部分支将两个长的外部分支分开，parsimony倾向于推断出不正确的树，如fig3的a-b，并且将长分支归为一组。 当T1中的分支长度足够极端时，支持正确树T1的位点模式xxyy的概率可能小于支持错误树T2的位点模式xyxy。因此， 序列中的位点越多，就越有可能在比xyxy更少的位点观察到xxyy模式，也就越有可能选择错误的T2树作为Maximum parsimony树 。 long-branch attraction已经在许多真实的和模拟的数据集中被证明了，这是由于parsimony不能改正在同一地点的多个变化或也不能适应两个长分支上的并行变化。

注意，如果假设的模型过于简单且忽略了位点率变化的话，基于模型的方法(即distance、likelihood和贝叶斯方法)也会遭受 long-branch attraction。 在重建深系统发育过程中， long-branch attraction(以及不同物种间核苷酸或氨基酸频率的不均等)是系统错误的重要来源。 在这种分析中，建议使用现实的替代模型和likelihood或贝叶斯方法 。 密集的分类单元样品，对于打破 long-branch attraction，并去除快速进化的蛋白质或位点也可能有帮助。

005.png

Maximum likelihood

M.1 Basis of maximum likelihood

最大似然方法是R. A . Fisher在20世纪20年代提出的一种估计模型中未知参数的统计方法。 它代表了数据中关于参数的所有信息。 参数的最大似然估计值(MLEs)是使可能性最大的参数值。 很多时候， 利用迭代优化算法，数值求解最小二乘问题。 MLEs具有理想的渐近(大样本)属性:它们是无偏的、一致的(它们接近真实值)和有效的(它们在无偏估计中方差最小) 。

M.2 Maximum likelihood tree reconstruction

由于计算能力和软件实现的提高，以及越来越现实的序列演化模型的发展，该方法得到了广泛的应用。 请注意，极大似然树估计涉及两个优化步骤:优化分支长度以计算每个候选树的树分数，以及在树空间中搜索最大似然树。 从统计的角度来看，树(拓扑)是模型而不是参数，而给定树的分支长度和替换参数是模型中的参数。 因此，极大似然树推理相当于比较许多统计模型，每个模型具有相同数量的参数。 上述MLEs的渐近性质适用于给出真树时的参数估计，不适用于极大似然树 。 给出真树时，MLEs具有吸引性的渐近性质，但不适用于极大似然树。

在距离计算时用到的所有替换模型都能在这使用。 用似然法对许多序列进行联合比较，使得适应更为复杂的序列演化模型成为可能。 分子系统发育学中使用的大多数模型都假设序列中各个位点的独立演化，因此似然是不同位点概率的乘积。 在任何特定位置的概率是祖先节点上未观察到的字符状态的平均值。 Likelihood和parsimony在这方面是相似的，尽管parsimony只使用最优的原始状态，而Likelihood对所有可能的状态取平均值。

早期的最大似然算法的实施包括PHYLIP，MOLPHY，PAUP4.0。现在的最大似然算法应用程序包括PhyML,RAxML,GARLI，不仅计算能快，而且在找高质量的似然分数的树上更有效。 最近在MEGA 5中加入了极大似然值，这使得没有计算机使用经验的生物学家更容易使用这种方法。

M.3 Strengths and weaknesses of the maximum likelihood method

最大似然法的一个优点是它所有的模型假设都是明确的，以便对它们进行评估和改进。 在Maximum likelihood 和贝叶斯方法中丰富的复杂进化模型的可用性是其相对于maximum parsimony的主要优势之一。 利用保守蛋白的现代推论几乎完全依赖于Maximum likelihood和贝叶斯方法。

如果目的是了解序列演化的过程，极大似然法比距离法或maximum parsimony有明显的优势。 likelihood ratio test可用于检验进化模型的拟合性，并检验有趣的生物学假设，如影响蛋白质进化的分子钟和达尔文选择。

最大似然的主要缺点是likelihood计算，特别是似然准则下的树搜索，计算量大。 另一个缺点是，如果模型是错误的，那么该方法可能具有较差的统计特性，这对贝叶斯分析也适用。

Bayesian methods

B.1 Basis of Bayesian inference

贝叶斯推理是一种通用的统计推理方法。 与极大似然不同的是，模型中的参数被认为是具有统计分布的随机变量，而在极大似然中它们是未知的固定常数。 在对数据进行分析之前，为参数分配一个先验分布，并与数据(或likelihood结合产生后验分布。 所有关于参数的推论都是基于后验分布的。 在过去的二十年中，由于计算方法的进步，特别是Markov chain Monte Carlo algorithms(MCMC算法)的进步，贝叶斯推理得到了广泛的应用。

B.2 Bayesian phylogenetics

开发了更高效的MCMC算法，消除了时钟假设(允许在无根树上有独立的分支长度)，并发布了MrBayes程序，这使得这种方法在分子系统学家中很流行。

006.png 007.png

B.3 Strengths and weaknesses of the Bayesian inference method

似然方法和贝叶斯方法都使用似然函数，因此具有一致性和有效性等统计特性。然而， 极大似然推理和贝叶斯推理是统计推理的对立哲学。 因此，贝叶斯推理的相同特征可以被看作是优点，也可以被看作是缺点，这取决于一个人的哲学。

首先，贝叶斯统计被认为是直接回答了生物学问题，并且产生了容易解释的结果：给定数据和模型后，树的后验概率就是树的正确的概率。 像似然分析中的置信区间这样的概念有一种人为的解释，许多统计数据的使用者无法理解。 在系统发育学中，定义树的置信区间是不可能的。 根据真实数据集计算的树和枝的后验概率往往过高。

其次， 先验概率允许包含关于树或参数的先验信息。然而，这类信息很少可用，而对先验的规范通常是一个负担用户;几乎所有的数据分析都是使用计算机程序中的“默认”先验进行的。

008.png

Statistical assessments of phylogenetic methods

系统发育推断的目的是估计树的拓扑结构和可能的分支长度。使用四条准则来判断树构建的方法。

1.Consistency. 当数据量趋近于无穷时，如果估计值聚集于真实参数值，则称估计方法是一致的。 当序列中位点数目增加时，估计值聚集于真实值时，则树重建方法是一致的。 如果假设的模型是正确的，那么基于模型的方法(即distance matrix、maximum likelihood和贝叶斯推理)是一致的。在某些模式下，Parsimony可能是不一致的;费尔森斯坦对此的论证引起了激烈的讨论。

2.Efficiency. 在参数的统计估计中，方差较小的无偏估计比方差较大的无偏估计更有效。 在系统发育学中，可以通过在给定位点数量的情况下恢复正确树或子树的概率来衡量效率。 树重构的复杂性意味着MLEs的渐近理论不再适用。 计算机模拟通常发现，在恢复正确的树时，maximum likelihood的efficiency比maximum parsimony或neighbour joining更高。

3.Robustness. 如果一个方法给出了正确的答案，即使它的假设被违背了，它仍然是具有稳健性的。 随着测序数据的快速积累，极大地降低了树重构中的样本误差，因此该方法的系统误差或稳健性变得越来越重要。

4.Computational speed. 这个属性很容易评估。Neighbour joining使用一种聚类算法来画树，并且非常快。 在maximum evolution、maximum parsimony和maximum likelihood等条件下寻找最佳树的方法比较慢。 贝叶斯方法的计算速度取决于链的长度(由MCMC算法生成)，这是高度依赖于数据的。 由于系统进化likelihood的计算成本较高，因此maximum likelihood和贝叶斯推理的速度通常要慢于maximum parsimony。 计算算法的巨大进步使得基于likelihood的方法对大数据集的分析成为可能。

Phylogenomic analysis of large data sets

P.1 Supertree and supermatrix approaches

对数百或数千个基因或蛋白质进行系统发育分析有两种方法，特别是当某些物种的某些位点缺失时。 supertree方法分别分析每个基因，然后使用启发式算法将每个基因的“子树”组合成所有物种的“超级树”。 分离分析有助于研究重建“子树”的差异或水平基因转移的流行程度。 然而，以所有基因为基础估计一个共同的系统发育是低效的。

在超矩阵方法中，将多个基因的序列串联起来，生成一个数据超矩阵，其中缺失的数据用问号代替，然后用这个超级矩阵重建树. 大多数超级矩阵分析忽略了基因间进化动态的差异。请注意，假设不同进化模型、不同的树和分支长度的超矩阵分析等同于分离或“超级树”分析。 当一棵共同的树是所有基因的基础时，理想的方法应该是对所有基因进行联合分析(超矩阵)，利用likelihood来适应进化过程中基因间的异质性。

P.2 Impact of missing data

许多基因组数据集是高度不完整的，因此大多数细胞在物种的基因矩阵将是空的。 虽然在理论上，似然函数(在极大似然和贝叶斯方法中)可以很好地适应缺失数据，但是这种大规模缺失数据和对齐间隙的影响还没有被很好地理解。模拟结果表明，在处理缺失数据时，最大似然推理和贝叶斯推理通常比neighbour joining或maximum parsimony表现得更好，贝叶斯推理表现得最好。 如果考虑到一些极端情况，即在去除排列gaps后，从不同的基因或位点计算出两两之间的距离，其中一些是快速进化的，而另一些则是缓慢进化的，则可以理解neighbour joining的糟糕表现。

P.3 Importance of systematic errors

在分析非常大的数据集时，几乎所有的bootstrap支持值或贝叶斯后验概率都被计算为100%，即使推断的系统发育可能在不同的基因间相互冲突，或者可能取决于使用的方法和模型。 因此，在这类分析中，系统误差比随机抽样误差要重要得多，而对违反模型假设的方法，即使它们的效率较低，也应该是可取的。（是否是放大系统误差，不得而知）

P.4 Data-partitioning strategies

数据划分的基本原理是将具有相似进化特征的基因或位点分组到相同的划分中，以便使用相同的模型描述同一划分中的所有位点，不同的划分使用不同的模型。 划分太细会增加计算时间，导致过度拟合，而划分太粗则会导致拟合不足或模型冲突。 然而，情况是复杂的，因为一些模型允许在取代率、氨基酸频率或取代模式的位点之间随机变化。 这种混合模型使用统计分布来适应不同位点间的异质性，而不需要进行数据划分。 选择使用分区模型还是混合模型是一个哲学问题:它分别对应于统计中对固定效应模型或随机效应模型的偏好。

目前的数据分割策略包括根据基因的相对替代率对基因进行分割，以及将编码基因的三个密码子位置分割成不同的分区。 likelihood比率测试也被用来决定两个基因是否应该在相同或不同的分区。总的来说， 数据划分与其说是一门科学，不如说是一门艺术，它应该依赖于我们对生物系统的知识:例如，假设所有基因都有相同的系统发育过程是合理的。