所用数据为一个属内不同种不同群体的叶绿体基因组序列，数量为80条。
发现用全长序列建树的时候，不适合选用太多外类群，否则ML法中会导致属内分枝的枝长特别短。原因应该是基因间隔区和内含子区域序列位点的差异较大。
枝长含义
NJ：表示遗传距离；
MP：性状状态变换的替换数；
ML/BI：该分枝上的相对进化数量（遗传变异量）；每个位点上的替换数（一般以每位点多少次核苷酸替换或氨基酸取代来表示）。
遗传距离
大多数情况以序列来说遗传距离就是两个OTU（个体、群体、物种或基因家族）之间序列的差异值。
序列比对
多序列比对用mafft得到的结果较为准确，muscle比对的速度较快。
多序列比对的绝大多数算法都是基于渐进比对的概念。简单来说就是先从两个序列的比对开始，逐渐添加新序列，直到所有的序列都加入为止。但是不同的添加顺序会产生不同的比对结果。所以由最相似的两个序列开始比对，由近到远逐步完成最为可靠。
mafft --thread 15 --auto 80-AcoeOut.fasta > 80-AcoeOut_aln.fasta
##比对时如果不清楚什么参数合适，加个参数 --auto，软件可以自动帮你处理
挑选保守位点进行下一步建树
序列比对完后，用于建树的序列位点必须保证具有良好的同源性。所以需要删除序列分歧很大的区域和gap区域。
我用的软件为Gblocks，主要目的是把有gap的位点全部去除，参数为-b5=n，其余的选项有-b5=h，h表示half 指去除在大于50%的序列中出现gap的位点。
Gblocks 80-AcoeOut_aln.fasta -t=d -b5=n
最大简约法（软件PAUP）
最大简约法的树长指所有性状在一棵树上的进化改变总数。
计算得到的结果可能会有许多树长相等的简约树，此时需要计算它们的一致树。分为strict consensus和semistrict consensus等，strict表示100%，在所有简约树中都出现的分枝，才会出现在一致树中，否则为梳子。这个阈值可以调。
一般文章中所用的系统树的拓扑结构都为ML或BI树，所以要把MP的bootstrap值标到ML/BI法的底树上。
1.cstatus ##查看数据特征，如简约信息位点数量等
2.tstatus
3.define outgroup
4.analysis-branch and bound ##数据量大用heuristic search
5.describetrees
6.contree ##计算一致树
7.bootstrap ##自举值
8.print bootstrap consensus
进化模型
DNA序列进化就是序列位点上的核苷酸随时间的变化，主要包括碱基替换、缺失和插入。
两条比对好的DNA序列的同源位点之间很容易看出碱基的相同或不同，但是在漫长的进化过程中实际发生了什么我们并不知道。最常见的当然是单次替换，但是当进化时间较长时，已经发生过替换的某些位点可能会再次发生替换，即多重替换。
DNA序列的进化模型将DNA的进化作为一系列随机突变来描述，并明确定义了4种碱基之间相互的替换速率。
DNA进化模型的参数主要有4类：
- 碱基组成频率
- 替换速率矩阵：指定了4种碱基之间在单位时间内或给定分枝长度下相互替换的概率。
- 不变位点比例
- 位点之间速率的变异：不同位点之间替换速率的差异。
ML法和BI法都需要选择合适的进化模型。模型选择软件具有的模型越多，检测结果越准确，但建树软件不一定支持该模型。
判断模型与数据拟合好坏的标准主要有AIC和BIC等。
最大似然法（软件IQ-TREE）
似然值是当模型（树和进化参数）为真时能够得到实际观测数据的概率。似然值是观测数据（即序列）的条件概率，其条件为计算似然值时依据的模型，而不是模型为真时的概率。
ML法建树的过程是先选择一个适合数据集的进化模型，然后对指定拓扑结构的一棵树优化分枝长度，以使得该拓扑结构的似然值最大化。通过计算不同拓扑结构树的似然值，将具有最大似然值的树看成是指定模型下的能够产生观测数据的最佳估计。
ML法采用的搜索方法主要是启发式搜索，步骤如下：
- 通过NJ树或逐步添加序列的方法构建初始树；
- 以初始树为基础通过各种分枝交换方法（TBR、SPR等）计算似然值，将最大似然值的树保存，并作为下一轮重排的初始树；
- 重复进行分枝交换，直到不能增加似然值为止。重排的最后获得的最大似然值树即为ML树。
建ML树的软件用RAxML的较多，但近来IQ-TREE的引用量一路上升。综合使用下来，个人感觉IQ-TREE的速度真快。
使用过程是下载了PhyloSuite的组件，从选模型到构树一站式操作还挺方便的。注意下载好后首先要配置用于不同分析的插件。
贝叶斯推论法（软件MrBayes）
BI法与ML法不同的是，前者根据提供的数据和选择的替代模型寻找可能性最大的树，而ML法则是寻找合适的树以使得数据的可能性最大。
推断系统发育树的步骤为：
- 选择一些树作为起始点；
- 判定这些树的似然值；
- 修改树的拓扑结构和分支长度；
- 计算出新树的似然值；
- 新树的似然值比旧树大，则接受新树。
如此就构成了一代，一次又一次的重复迭代，直到新树的似然值不再有明显变化，即树的似然值不再有显著区别，参数已收敛为止。如果没有收敛，适当的增加代树继续跑。
如何判断参数是否已收敛
软件运行完毕后，看结果文件的分离频率平均标准差值（Average standard deviation of split frequencies）
该值<0.01时，说明两次运行的结果差异很少，参数已收敛；
该值>0.05时，需要继续运行。
同样是在PhyloSuite中运行
网友评论