生物信息学初识篇——第二章：序列比对（2）

生物信息学初识篇——第二章：序列比对（1）

四、双序列局部比对

一长一短的两条序列，比较局部比比较全长更有意义。局部比对的算法和全局比对很相似，只是在选最大值时通过增加了第四个元素“0”，来达到比对局部的效果。序列p和序列q，一长一短，其他输入值跟全局比对的一样（图2.26）

图2.26 局部比对输入值：序列p和序列 q替换记分矩阵，空位罚分

局部比对的计算公式在全局比对的基础上增加了第四个元素“0”。得分矩阵初始值仍是0，但第一行和第一列与全局比对不同，全是0（图2.27）。

图2.27 局部比对计算公式及得分矩阵

从 s(1,1)开始要选择四个值中的最大值。除了上面格 s(0,1)+gap=0+-5=-5，左边格 s(1,0)+gap=0+-5=-5，斜上格 s(0,0)+w(1,1)=0+-3=-3，还有一个 0。max(-5,-5,-3,0)=0。并且这个 0 既不是从上面格，也不是从左边格，以及斜上格三个方向来的，而是来自于公式里增加的“0”，所以不用画箭头（图2.28）。

s(1,4)的计算：上面格 s(0,4)+gap=0+-5=-5，左边格 s(1,3)+gap=4+-5=-1，斜上格 s(0,3)+w(1,4)=0+0=0，还有一个0。max(-5,-1,0,0)=0。这个0和s(1,1)的0是不一样的。这个0应该画上斜上的箭头（图2.29），因为它可以来自公式中的 0，也可以来自斜上格。而s(1,1)的0没有箭头因为它只来自公式中的0。两种情况虽然都是0，但来源不同，一定要通过箭头标识清楚。

图2.28 得分矩阵中的 s(1,1)

图2.29 得分矩阵中的 s(1,4)

按照公式，填充满整个得分矩阵（图2.30）。与全局比对不同，局部比对的得分不是在右下角，而是在整个矩阵中找最大值。这个最大值才是局部比对的最终得分，他可能出现在任何一个位置。这次箭头追溯也不是从右下角到左上角，而是从刚刚找到的最大值开始追溯到没有箭头为止。追溯箭头终止的位置也可以是得分矩阵中的任何一个位置。

图2.30填满分值和箭头的得分矩阵

最后根据标记好的箭头写出比对结果（图2.31）。从左上到右下标记的红色箭头依次是：斜箭头字母对字母，C 对 C；斜箭头字母对字母，G 对 G。相比这两条序列的全局比对结果，两端的空位在局部比对中就全部被忽略掉了。

图2.31 序列 p 和序列 q 的局部比对与全局比对的比较

五、一致度和相似度

两条长度不同的序列做全局比对，然后计算全局比对中一致字符的个数和相似字符的个数，再除以全局比对的长度，就可以得到它们的一致度和相似度了。比如下面这两条序列：

首先做出它们的全局比对，比对中一致字符的个数是 4 个，全局比对长度 6，一致度=67%。相似字符个数 1，相似度就是(4+1)/6=83%。

把长度相同的两个序列计算一致度和相似度的方法重新规范一下。尽管长度相同，但是做出的全局比对的长度并不一定等于序列的长度，比如下面这两条序列：

上下各加入一个空位，全局比对的长度就不等于序列的长度了。所以不管两条序列长度是否相同，都要先对它们做全局比对。让两条序列先以最优的方式比对起来，再从全局比对中数出一致字符和相似字符的个数，除以全局比对的长度，来得到它们的一致度和相似度。

六、在线双序列比对

（一）、EMBL全局双序列比对工具

目前，使用率最高的是 EMBL 网站的双序列比对工具（http://www.ebi.ac.uk/Tools/psa）。打开页面，上面有全局比对工具、局部比对工具、还有基因组比对工具。首先看全局比对中的蛋白质序列比对工具（图2.31）。输入值非常简单，把要比较的两条蛋白质序列贴在输入框里或者上传。如果想要进一步设置比对的参数，可以点 More options。从这里可以选择使用哪种替换记分矩阵。按照之前讲过的原则，选择 PAM 矩阵或 BLOSUM 矩阵。如果实在不知道选哪个矩阵，就闭着眼睛选 BLOSUME62，下拉菜单里默认选的就是BLOSUM62。除了选择替换记分矩阵，这里还可以设置空位罚分，也就是gap的分值。这里实际上是让你选空位对字母的情况罚几分，所以显示的是正数，但在计算的过程中还是按照负数处理。这里的 gap 分好几种，一种叫“gap 开头（GAP OPEN）”，另一种叫“gap延长（GAP EXTEND）”，gap 开头就是连续的一串 gap 里面打头的那一个。gap 延长就是剩下的那些gap，这一串里，第一个gap 是gap开头，后面的都是 gap 延长。单独的一个gap按gap开头算。gap开头和gap延长可以分别定义不同的罚分。默认情况下，gap 开头罚分多，gap 延长罚分少（图2.32）。当gap开头小，gap延长大的时候，做出来的比对里面，gap很分散，极少有连续长串的gap出现（图2.33-A）。开头的一串gap是个例外，因为 seq2 太短，seq1的这一段只能跟gap相对。其他部分的gap都是分散出现的。这和我们默认参数（gap开头大，gap延长小）做出来的比对结果是截然不同的（图2.33-B）。在实际应用中，需要根据不同的情况选取不同的 gap 罚分，以满足不同的生物学意义。如果你对结果没有什么预期，那就请保持默认的参数。

除此之外，结尾的 gap 也可以划分出不同的种类并赋予不同的罚分，如果把 END GAPPENALTY 选成 true，就可以设置结尾的 gap 罚分了。结尾 gap 不太常用，特别是在做亲缘关系较近的序列比对时，是否设置结尾 gap，比对结果差别不大。

图2.31 EMBL全局双序列比对输入页面

图2.32 EMBL全局双序列比对输入页面

图2.33 设置不同的 GAP OPEN 和 GAP EXTEND

（二）、EMBL局部双序列比对工具

EMBL的局部双序列比对工具可以选择经典的 Smith-Waterman 算法。More options 里面的参数设置和全局比对是一样的。在这个例子里，我们保持所有参数都为默认值，点提交（图2.34）。

图2.34 EMBL 局部双序列比对输入页面

从比对结果可以看出（图2.35），只有中间黑色的相似的部分出现在比对结果中了，两头红色的不相似的部分被忽略掉了。也就是只返回了局部最相似，得分最高的片段的比对结果。

如果给这两条序列做全局比对的话，会发现，绝大部分位置对得都很差，只有中间这一段对的还不错（图2.36）。所以，有时候两条序列并不同源，它们只是有一个功能相似的区域，这时用局部比对我们就能很快找到这一区域在两条序列中的位置。但是如果做全局比对的话，结果就不如局部比对明显了。

图2.35 局部双序列比对结果

图2.36 全局比对与局部比对的比较

（三）、其他在线双序列比对工具

可以做双序列比对的工具很多（图2.37）。不同网站都有自己的比对工具，所使用的算法也不尽相同，但是它们的核心算法都是讲过的 Nidelmann-Wunsch 和 Smith-Waterman 算法，只是在他们的基础上有所变化，有所升级。

Biotools 的双序列比对工具无论是核酸序列还是蛋白质序列都能做。全局比对，局部比对，还包括蛋白质二级结构辅助的序列比对。功能比 EMBL 的只多不少。当然这都不是最关键的。关键的大招是除了能给出序列比对，还能给出得分矩阵。

图2.37 其他在线双序列比对工具