提到构建系统发育树，我们大概都会想到一堆的软件，如下表：

(图片来源：H.O. Minoarivelo, 2015)

这些软件背后的算法主要分为四种：比如当下最广泛使用的MEGA，就提供了最大似然法、最大简约法、距离矩阵法，可以选择bootstrap法测验，每个方法也提供了不同的进化模型以供选择。而近年来引用量逐年攀升的MrBayes，采用的则是贝叶斯推断法。

这四种方法到底有什么不同，我们来看一看。

首先是距离矩阵法(distance matrix method)

顾名思义，是基于距离矩阵构建的系统发育树，那么什么是距离矩阵？

wiki给出的解释是a distance matrix is a square matrix (two-dimensional array) containing the distances, taken pairwise, between the elements of a set. 即距离矩阵是一个二维的数组，这个二维数组中包含了数据集中两两元素之间距离。

但是！ 距离是怎么计算的呢？

在系统发育分析中，我们称为遗传距离，指的是所有成对实用分类单元OTU(operational taxonomic units)之间的距离(或称差异值)，实用分类单元可以是个体、群体或物种、种上分类单元或基因家族。所选用的数据可以是序列也可以是其他性状数据。

如下，假定有n个OTU，那么将组成n x n 的矩阵，需要计算的成对距离有n(n-1)/2个。

假设，我们分析的数据是序列数据，要找出OTU之间的相似性或者差异性，那么我们自然会先对它们进行比对，比对之后，我们发现，有单个碱基变异的位点，也有插入缺失的碱基，那么我们应该选用什么样的值来表示它们的差异性或者遗传距离呢？

• 总替换分歧度：根据序列之间所有类型的差异计算的未校正的遗传距离

n为比对后的序列长度，nd为两条序列之间不匹配的位点数目。

也可以从相似系数中计算，首先计算比对序列之间的相似性值S，然后将S转化成差异值D(即距离)

ni和M为序列比对后具有相同残基的位点数；n为所比较序列的位点总数；L为所比较序列的有效长度；U为非匹配位点数(不含空位)；G为所引入空位的位点数；W为对空位的加权值(0,1)，一般取中间值0.5。

(相同残基位点，严格意义上是指具有相同的碱基或者氨基酸位点)

• 转换/颠换分歧度：序列对之间发生的转换位点数/序列长度，或者颠换位点数/序列长度；

• 同义/非同义替换分歧度：同义替换距离ds；非同义替换距离dn，适用于计算蛋白的分歧度。

• indel距离：空位距离

• 其他的一些分子数据距离

每个距离的计算都有涉及到一些模型，其中，单个碱基变异的模型研究的比较清楚，目前主要用SNP(单碱基变异多态性)数据构建系统树。

Jukes-Cantor单参数进化模型和Kimura两参数进化模型

模型的构建基于假设，在分子进化研究中，往往假设序列是同源的，它们具有单一的祖先序列，且这一祖先序列在进化过程中发生了一系列核苷酸的突变。

在该假设前提下，Jukes和Cantor(1969)进一步假设每种碱基具有同等概率突变为另三种碱基，频率常数为μ/3，μ为碱基替换频率。

而kimura(1980)考虑到转换 (同类型碱基的替换)和颠换(不同类型，如嘧啶替换为嘌呤) 具有不同的发生概率(α和β)，提出了两参数模型。

距离法主要包括3种主要方法：

1：非加权平均连接聚类法或UPGMA法

UPGMA是unweighted pair-group method using an arithmetic average,

翻译过来就是应用算术平均数的非加权成组配对法，这种方法虽然名字很长，但是主要内容就是将成对距离的平均值作为两两成员之间的距离，是早期的一种聚类方法。

举例说明，给出人类、黑猩猩、大猩猩、猩猩、长臂猿的线粒体序列，比对统计两两之间的核苷酸替换数，采用JC69模型估计遗传距离。

(从左到右依次是什么哈哈哈哈，猜一猜)

已知，人类和黑猩猩的距离最近，那么将人类和黑猩猩聚为一个新的OTU，大猩猩与新OTU的距离：1/2(人类和大猩猩的距离+黑猩猩和大猩猩的距离)，依次可以计算出其他类别和新OTU的距离，得到一个新的距离矩阵。

此时分类单元变成了4个，其中又继续将距离最近的两个聚为一个新的OTU，重复上面的过程，直至最后剩下两个OTU为止。

最后，将聚类的结果以树状图的形式呈现，两个分支的长度之和表示序列的遗传距离。

2：Finth-Margoliash法

UPGMA法的假设前提是所有分支的突变率是相等的，为一个固定的常数。

跟UPGMA将两个距离最近的类别聚为一个新OTU类似，F-M法对距离最近的两个类别，引入节点作为最近共同祖先，比如人类A，黑猩猩B，其他归为X，引入AB最近共同祖先C，求解三个二元一次方程，得出A，B到C的距离，以及节点C的高度。依次类推，最终求解出所有所有类别到节点的长度，得出无根树。此时，不必要求分支具有相同的突变率。

但如果需要得出有根树，需要依据树根到G和H的距离相等，才能得出树根的位置，此时，突变率已经默认相同了。

3：邻接法(NJ法)

邻接法应该是这三种方法中，我们最熟悉的一种了。

邻接法是Saitou和Nei(1987) 根据最小进化原理创建的一种叠加树的方法，和上面两种方法相似，邻接法不断地将距离最近(相邻的点)合并成新的点，来使树的总距离最小。其不同的地方在于，对距离进行了校正。原始距离用于估算距离矩阵，校正距离用于判断节点之间的连接顺序。

1：计算每个末端分类单元的净分歧度。

比如，第i个的净分歧度

简单直白的说，就是将人类到黑猩猩，大猩猩，猩猩，长臂猿的距离加起来作为人类的净分歧度。

2：计算最小速率校正距离。

前面我们提到过未校正的距离，比如上面根据JC69和K80模型计算的距离就是未校正的距离。由于分歧时间不断延长，某些快速进化位点很可能发生了多重替换，而多重替换后将消除之前的进化痕迹。此时，就需要对原始距离进行校正，而最关键的就是需要一个涉及位点固定的速率矩阵。

3：将校正距离最小的两个节点i，j合并为新节点u，计算节点i，j到u的距离

如此，获得新节点u，以及剩下的节点，节点数N-1，继续三步走，直到只剩下两个节点，系统树生成。

小结一下：

距离法是以距离矩阵为基础的分析方法，对原始数据进行遗传距离估计，此时的距离也称为观测距离，观测距离在大多数情况下是被低估的，所以要根据进化模型进行校正，得到校正距离。得到距离矩阵后，如果距离关系满足超度量关系，则可以直接使用聚类分析构建超度量树(比如UPGMA)，否则就要使用邻接法，最小进化法或叠加法等。

UPGMA要求符合分子钟的假设，即各个分支的突变速率相等，而大部分的数据可能并不符合这个假设，所以相比之下，NJ法不受分子钟的限制，在准确性和运算速度上更好一些。普遍来说，距离法已经不常用了，但是可以将NJ法可以作为ML树的初始树，以便节约时间确定模型，参数等。

补充：

什么是超度量关系：

生成系统树之后，我们仍需要判定它是否可靠，所以需要对其进行可靠性检验，比如常用的自举法检验(bootstrap)，什么是自举法检验，介绍完几种建树的方法之后再写一写吧。

详细请参考：

黄原 《分子系统发生学》 第八章，2012

樊龙江 《生物信息学》 第1-6章，2017

https://en.jinzhao.wiki/wiki/Distance_matrix