比较不同长度序列的方法(打点法与动态规划算法)
紧随上文,知道了如何比较两个长度相同的序列,那么如何比较不同长度的序列,并求出它们的一致度与相似度呢?
一般有两种方法:1)打点法;2)序列比对法。
打点法
打点法非常简单,理论上在纸上都可以进行。
实际上就是将两个序列间的点两两组合。
若A序列长10,B序列长9。
则需要进行10x9 = 90 次的比较
解读打点法信息
打点法序列的对角线及对角线的平行线,代表两条序列中相同的区域。
THEFA; TCAT; AT
-
比如下面这两个序列,可以说是风马牛不相及。
给自己打点
- 除了可以用打点法给不同的序列打点,还可以用一条序列给自己打点。从而可以发现序列中的重复片段。
这样的打点矩阵必然是对称的,并且有一条主对角线。
在横向或者纵向上,与主对角线平行的短平行线所对应的序列片段就是重复的部分。
- 还可以通过这种方式,发现串联重复序列。(tandem repeat)
短串联重复序列,short tandem repeat, STR, 也叫做微卫星DNA,是一类广泛存在于真核生物基因组中的DNA串联重复序列。由2-6bp 的核心序列构成,重复次数通常在15-30次。STR 具有高度多态性,即存在重复次数的个体间差异,而且这种差异在基因遗传过程中一般遵循孟德尔共显性遗传规律,所以被广泛应用于法医学个体识别、亲子鉴定等。
在线打点软件
dotlet 软件
提供序列录入,以及窗口大小和替换记分矩阵的选择。
-
因为使用的edge 浏览器不支持java,可好在dotlet 还提供了JS版本。
开始打点
这里可以使用dotlet 提供的seq1及seq2 序列文件。
并调整相关参数。
window size 与scoring matrix 分别对应以多长的序列为单元打一个点,以及替换记分矩阵的类型。
若window size 为1, 则代表一次只比较1个字母(同之前打点举例)。若为15(默认),则代表一次比较15个字母,通过比较它们整体的相似度来确定是否打点。
需要注意的是,这种比较并不是1-15,16-30 这样的方式来进行比较。而是1-15,2-16,这样的方式。
因此选择 1还是n(比较n个字母),所比较的次数,也只是相差了 n-1。
这里需要注意,一般不同序列的比对,使用BLOSUM 62。或者也可以根据序列的差异程度,自行选择其他BLOSUM 或PAM矩阵。
只有当比较两列相同序列时,才选择identity 矩阵类型。
解读打点结果
将同一序列进行自己打点。
由图可见,和手工画图一样,也是在白纸上打下黑点。
而所有自己和自己打点的矩阵都有的共性一样,也有一条分隔两边的主对角线。
使用鼠标点击打点图,便可以定位在该位置。通过下方的序列信息显示区域,相同的序列信息也会高亮出来。
而序列信息显示区域红色的方框,即是之前选择的window size的数值。默认为15。红色方框中间的字母,也上下延伸出定位信息,表明它们在各自序列上的位置。
在java 版本中dotlet 还提供丰度信息。
分值由-60-165,分值越高,相似度越高,为x轴;纵轴为得该分值得点的个数。
可以通过调整灰度条,屏蔽相似度低的点。突出对比。
若是自比较打点,包括主对角线在内,一侧的平行对角线有几条,序列中就有几条重复序列。大的重复序列中也可能包含小的重复序列。
序列比对法
使用打点法dotlet,只能让我们大致了解两条序列的相似程度(对角线),但无法精确描述。
如果需要精确描述两条序列的相似程度,则需要使用序列比对法。
序列比对,aligment,也叫对位排列,联配,对齐等。是运用特定的算法找出两个或多个序列之间产生最大相似度得分的空格插入和序列排列方案。是序列比较的标准方法。
序列s和序列t的比对:
其主要方式就是通过将序列上下排列在一起,并通过在某些位置插入空格(gap),比较上下字符的匹配情况,从而找出使两条序列产生最大相似度得分的排列方式和空格插入方法。
序列比对可以分为:1)多序列比对;2)双序列比对。
而双序列比对又可以分为全局比对与局部比对。
局部比较比的是对的好的局部,对的不好的部分将会被忽略不计。
双序列比对
双序列比对的全局比对及其算法。
现在所有的序列比对的算法,几乎都是从Needleman-Wunsch 算法衍生而来。
输入值与替换记分矩阵
序列p:ACGTC
序列q:AATC
相比相同长度序列比较中提到的BLAST 矩阵等,Needleman-Wunsch 算法中多出来了gap,且记为-5分。
而空位对空位,是毫无意义(嗯哼?没事找事?),完全可以剔除掉。
得分矩阵
得分矩阵由序列 p,q构成。
并在序列开始位置各留一个空行与一个空列。
记分矩阵的算法
开始记录
1)首先从第零行写起
因为q 在列上,而第零行中,q 均为0,所以按照
s(0, j) = gap * j
计算。这实际是一种极端情况的假设,即序列p 中的字母,全部对空位时的得分。
2)接下来填第零列
和第零行一样,按照
s(i, 0) = gap * i
计算,第零列的意义也是一种极端情况的假设,即序列q 中的字母,全部对空位的得分。
3)接下来填s(1, 1)
这时候的情况稍微复杂了一些。
需要填入的值有三种选择的情况。
第一,是s(0, 1) + gap
。其表示序列p(行上)对空位。得分为-10。
第二,是s(1, 0) + gap
。其表示序列q(列上)对空位。得分为-10。
第三,是s(0, 0) + w(1, 1)
。其表示序列q 与序列p 在相应位置发生匹配。得分为10。
因此选择三种情况得分的最大值。
4)用箭头记录轨迹
我们需要通过箭头,记录最终选择的数值来源的箭头。
5)小练习。
将下面的动态规划得分矩阵用箭头标注完。
~ 揭晓答案
6)确定比对序列
箭头绘制好了以后,右下角的得分也就是最终得分矩阵得分的最大值。将右下角箭头回溯至起始位置,并将这些箭头标记出来。
根据3)中的计算我们知道,向上的箭头,表示序列q(列)对应的序列p(行)为空值;而向左的箭头,表示序列p(行)对应的序列p(列)为空值;斜方向表示序列p 与序列q 相对。
因此,顺着最大记分的箭头,便可以画出最终的全局比对结果。
局部比对
上边的内容,为双序列比对中的全局比对(global alignment),是适用于两个长度近似的序列。
那么对于两个序列长度不等,且差异很大,该怎么办呢?
试试看接着用全局比对
小实验
假设有三个序列a,b,c
很明显,当使用全局比对时,a与c 的得分显著低于 b与c 的得分,那是不是可以说a与c 更相似呢?
如果我们单纯只看红框里的部分的话,不难发现,a 序列与c 序列是更相似的,而且有很高的一致度。通过计算,也验证了猜想,单纯比较序列a 红框里的部分碱基,与c 的得分更高,为24。
这种通过比较两个序列局部位置的比对方式,称为局部比对(local alignment),一般也被用来比较一长一短两个序列。
局部比对与全局比对的差别
从上面的解释也不难看出,局部比对可以说是全局比对的衍生物。因此它们的算法也很相似。只是在选择最大值时多了一个零,来达到局部效果。
开始比对
假如我们的手头上有两个序列
序列p:ACGTC
序列q:CG
和全局比对一样,我们开始画表。
只是这一次,多了0这个选择。
接着,按照之前的经验,我们可以把表格填满。
与全局比对不同,局部比对最终的得分不是在右下角。而是在整个矩阵中,找最大值。它可能出现在任何位置。
而箭头追溯也不是到左上角,而是从最大值开始,追溯到没有箭头为止,也可以是矩阵中的任意一个位置。
而空位也就全部被忽略不记了。
总结
因此无论序列是否长度相同,都需要做全局比对。