一、简介
1.1 介绍
MUSCLE(Multiple Protein Sequence Alignment)是一款蛋白质水平多序列比对的软件，在速度和精度上都优于 ClustalW。在进行多序列比对的时候，大多数情况下可以优先使用 MUSCLE。有本地版和在线版，在线版网址如下：http://www.ebi.ac.uk/Tools/msa/muscle/。
1.2 算法原理
算法：MUSCLE 先使用渐进式比对(progressive alignment)获得初始的多序列比对，再使用横向
精炼(horizontal refinement)迭代提高多序列比对结果。
1. 使用数串(k-mer counting)方法构造序列间的全局比对和局部相似度
2. 填充序列间距离的三角矩阵
3. 使用 UPGMA 或 NJ 法构建序列发生树，并确定无根树的根
4. 从叶节点开始向上推测父节点的渐进式比对，最后产生根节点的多序列比对
5. 根据得到的多序列比对，计算任两序列间的相似度
6. 计算 Kimura 距离矩阵，构建发生树
7. 比较新生成的树和原来树的差异，如果有节点的重排，跳转到步骤 4
8. 从树上砍断一个枝，产生两个子树，每次砍断的位置是按和根的距离降序排列的
9. 分别计算两个子树的多序列比对，并对两个结果比对得到新的多序列比对
10. 如果新的比对结果的 SP 分数(sum of pairs)降低，保留这个新的比对结果，反之丢弃。反复迭代 8－>9－>10，直到分值不再降低或达到最大迭代次数。
1.3 使用'muscle'做多序列比对
MUSCLE做多序列比对十分简单，大多数情况下用户只需要指定输入输出文件即可。
例如，输入文件'file1'为：
>a
ATGAGGTAGAGATAGCCGG
>b
ATGGTTAGCCGG
运行命令：
$ muscle -in file1 -out file2
则输出结果'file2'为：
>a
ATGAGGTAGAGATAGCCGG
>b
ATG———-GTTAGCCGG
运行程序log：
MUSCLE v3.8.31 by Robert C. Edgar
http://www.drive5.com/muscle
This software is donated to the public domain.
Please cite: Edgar, R.C. Nucleic Acids Res 32(5), 1792-97.
1 2 seqs, max length 19, avg length 15
00:00:00 11 MB(-1%) Iter 1 100.00% K-mer dist pass 1
00:00:00 11 MB(-1%) Iter 1 100.00% K-mer dist pass 2
00:00:00 12 MB(-1%) Iter 1 100.00% Align node
00:00:00 12 MB(-1%) Iter 1 100.00% Root alignment
1.4 生物信息中的利用
目前muscle主要用来在基因组进化部分，因为构建进化树和计算选择压力，都需要将序列对齐，muscle小而快，毋庸置疑是最好的选择。
二、构建进化树
注：本次操作使用的muscle版本为：MUSCLE v3.8.1551
muscle构建进化树一般需要两步，输入文件'test.fa'：
>1
CAGACAGTTTATAGCAACACTGCCAACTCAGCTCCCGCCACAGAAATCGGAGGGAGGGTT
CTGGAAAACACGAGTTCAAACTGGGAGCACCCCAGATGTGGG
>2
AGTTTATAGCAACACTACCAACTCAGCTCCCGCCACAGAAATCGGAGGGAGGGTTCTGGA
AAACACAGGTTCAAACTGGGAGCACCCCAGATGTGGG
>3
CAGGCAGTTTATAGTAACACTACCAACTCAGCTCCCGCCACAGGACCGGAGGGAGGATTC
TGGGAAAAACACGGATTCAAACTGGGAGCACCCCAGATGTGGG
>4
AGTTTATAGCAACACTACCAACTCAGCTCCCGCCACAGAAACCGGAGGGAGGGTTCTGGA
AAACACAGGTTCAAACTGGGAGCACCCCAGATGTGGG
>5
AGTTTATAGCAACACTACCAACTCAGCTCCCGCCACAGAAATCGGAGGGAGGGTTCTGGA
AAACACAGGTTCAAACTGGGAGCACCCCAGATGTGGG
2.1 多序列比对
运行命令：
$ muscle -in test.fa -out test-aligned.fa
输出结果'test-aligned.fa'为：
>3
CAGGCAGTTTATAGTAACACTACCAACTCAGCTCCCGCCACAG-GACCGGAGGGAGGATT
CTGGGAAAAACACGGATTCAAACTGGGAGCACCCCAGATGTGGG
>1
CAGACAGTTTATAGCAACACTGCCAACTCAGCTCCCGCCACAGAAATCGGAGGGAGGGTT
CTGG--AAAACACGAGTTCAAACTGGGAGCACCCCAGATGTGGG
>2
-----AGTTTATAGCAACACTACCAACTCAGCTCCCGCCACAGAAATCGGAGGGAGGGTT
CTGG--AAAACACAGGTTCAAACTGGGAGCACCCCAGATGTGGG
>5
-----AGTTTATAGCAACACTACCAACTCAGCTCCCGCCACAGAAATCGGAGGGAGGGTT
CTGG--AAAACACAGGTTCAAACTGGGAGCACCCCAGATGTGGG
>4
-----AGTTTATAGCAACACTACCAACTCAGCTCCCGCCACAGAAACCGGAGGGAGGGTT
CTGG--AAAACACAGGTTCAAACTGGGAGCACCCCAGATGTGGG
2.2 构建进化树
运行命令：
注：这里的输入是已经比对好的结果文件test-aligned.fa
$ muscle -maketree -in test-aligned.fa -out test.newick
这里默认是使用'UPGMA'方法来构建的，如果要使用'NJ'方法可以这样：
$ muscle -maketree -in test-aligned.fa -out test.newick -cluster neighborjoining
输出结果'test.newick'（'newick'格式）：
(
3:0.0280448
,
(
1:0.0160781
,
(
(
2:-0
,
5:-0
):0.00519213
,
4:0.00519213
):0.010886
):0.0119667
)
;
2.3 'newick'格式说明
'newick'格式由英国著名数学家`Arthur Cayley于1857年提出，目的是用数学的方法表示树形结构。标准的'newick'格式使用一对小括号表示一个'interior node'，可打印的字符（不能包含：空格、冒号、分号、小括号、中括号）表示节点，最后以一个分号结束。
例如：
(B,(A,C,E),D);
表示：
另外需要注意几点：
1. 节点名称可以为空，例如：(Alpha,Beta,Gamma,Delta,,Epsilon,,,);
2. 可以在节点名称后面加上分号和一个数字来表示分支长度，例如：(B:6.0,(A:5.0,C:3.0,E:4.0):5.0,D:11.0);
3. 除了节点名称和分支长度之外的任何地方都可以加入空白，例如：
(
B:6.0,
(
A:5.0,
C:3.0,
E:4.0
):5.0,
D:11.0
);
2.4 使用'biopython'画图
注：需要安装biopython和matplotlib，本次操作使用的版本为：biopython==1.72、matplotlib==3.0.0
有了muscle输出的newick格式的结果文件，我们就可以可视化了，可视化的方法非常多，这里使用biopython可视化：
from Bio import Phylo
tree = Phylo.read('test.newick', 'newick') # 如果有test.newick里面有多个树，则使用Phylo.parse方法读取
Phylo.draw(tree)
结果为：
更多生信分析需求请加微信：13120220117