仇瑛姿

摘 要：基因序列匹配是通过匹配算法，寻找两个或者多个序列之间的相似性和同源性，常见于对蛋白质序列和DNA序列的匹配。基因序列匹配算法从始至今，在保证匹配结果精度的前提下，不断演变出匹配时间较少的算法，以提供更优的匹配性能。文章即是介绍基因序列匹配算法的演变和实际的应用。

关键词：基因匹配；算法演变；算法应用

中图分类号：Q343.1 文献标志码：A 文章编号：2095-2945（2018）13-0063-02

Abstract： Gene sequence matching is to find the similarity and homology between two or more sequences by matching algorithm， which is usually used to match protein sequence and DNA sequence. From the beginning to now， gene sequence matching algorithms have evolved to provide better matching performance under the premise of ensuring the accuracy of matching results. This paper introduces the evolution and practical application of gene sequence matching algorithm.

Keywords： gene matching； algorithm evolution； algorithm application

绪论

随着生物科学和计算机科学的迅猛发展，计算机和生物相结合形成一门新学科，利用计算机对数据的快速处理能力，挖掘大量而复杂的生物數据，为生物医学的飞跃提供了有利条件。从1990年人类基因组计划至今，已完成了人类基因组DNA30亿个碱基对的测序，破译了人类超过95%的遗传信息，我国也在2017年12月启动了“中国十万人基因组计划”，通过绘制中国人精细的基因组图谱，来研究疾病健康和基因遗传的关系。世界各国对生物基因组测序工作极快地展开，而且可以预测，今后DNA测序需求的增长将更为惊人，这些海量数据的积累和运算，挖掘与分析，离不开现代的计算机技术。

由于生物基因库的不断扩充和细化，提高基因序列匹配的性能迫在眉睫。基因序列匹配算法也由原来单一的动态规划，全局匹配逐渐发展成启发式的局部最优算法。同时为了获得更准确的数据，在查询匹配序列之前，添加了对基因序列的过滤条件，算法的不断创新加之对前人算法的改良，使基因序列匹配算法逐步走向成熟。

1 全局匹配算法

全局匹配即整体考虑整条基因序列相似性，一条序列转化成另一条序列需要的最小距离决定了该序列相似性的得分值，距离越小，则得分越高，相似度也越高。较为常见的是动态规划算法，该算法通过拆解问题，之后处理问题之间不同状态之间的关系，使问题可以通过递归的方式去分解[1]。

动态规划算法适用于整体问题最优解包含子问题最优解和某一状态不会被前一状态影响，只与当前状态有关的场景。在基因序列匹配中常见的动态规划算法是：隐Markov模型（HMM）+Viterbi算法。

以下举例利用动态规划算法寻找序列相似性问题的方法（两条DNA序列的比较）[2]：

Step1：将两条DNA序列（Seq1， Seq2）的所有碱基对生成字符串矩阵，并得出两个序列所有碱基对应的对应关系的距离值（Seq1转换成Seq2的得分值）。计算两个碱基距离一般有两种方法：一种是汉明距离，即不考虑空位，移位的情况，一一对应；另一种是编辑距离，即将Seq1的一个位置上的碱基变换成Seq2的一个位置上删除/替换/插入的最少基本操作数目。

Step2：通过矩阵计算出从初始位置（两个序列的起点）走到矩阵下一个位置的距离（或者得分）。

Step3：遍历整个矩阵的位置后得到所有位置的得分。

Step4：从终点位置向前找距离最小的位置（或者得分最高的位置）。

隐Markov模型是一种重要的统计方法，在生物信息学中应用于线性序列分析或基因发现等发main，利用它对普通动态规划算法进行简化，这是由于普通的动态规划算法需要遍历矩阵所有位置，对于大量的基因匹配，会耗费相当长的时间，而HMM模型，包含隐含在两个状态之间的转移概率，可以减少遍历节点，加速普通动态规划的时间。

Viterbi算法应用到基因序列匹配[3][4]：

从HMM模型中可以推导出诸如核苷酸序列ACAC出现的概率：

P（ACAC）=P（A1）*P（2|1）*P（C2）*P（3|2）*P（A3）*P（4|3）*P（C4）

Viterbi算法要解决的就是HMM模型中从前t个碱基状态到最终的碱基状态k的观察结果中，取最后可能对应状态的序列，获得Viterbi路径。

2 局部匹配算法

局部匹配算法的核心是利用启发式在局部获得最优解，然后通过局部最优解得到全局最优解。这种方法精度没有全局动态规划的方法高，但比之前的方法快100倍左右，牺牲了少量精度换取了执行性能。要注意的是，动态规划算法其实也可以通过分治法转变成局部动态规划匹配，但这种算法不在本文讨论范围之内。常见的局部匹配算法有FASTA算法和BLAST算法。

FASTA算法和BLAST算法均用到了启发式，以简化运算步骤。启发式即利用先验概率获得预知能力，对不确定的状态进行判断，在某些极端情况下，启发式会获得很差的结果，但通常可以在合理的时间内获得正确的答案。启发式算法的结构是邻域搜索，即从一组初始解中产生若干邻域解，从邻域解中获得并更新当前的状态，然后迭代并取满足收敛条件的状态。

FASTA算法[5]是由Lipman和Pearson推導实现的，具体实现步骤如下：

Step1：从两个序列中寻找到长度为k的共同片段（k-tuple/ktup），诸如：

Seq1：GCATGCGCCCAC

Seq2：GAATGCGCCAAC

Seq1和Seq2的长度N1=N2=12，取k=2，则每一条序列都可以得到N-k+1个words，如下：

Seq1-words：GC（1）|CA（2）|AT（3）|TG（4）|GC（5）|CG（6）|GC（7）|CC（8）|CC（9）|CA（10）|AC（11）

Seq2-words：GA（1）|AA（2）|AT（3）|TG（4）|GC（5）|CG（6）|GC（7）|CC（8）|CA（10）|AA（11）|AC（12）

Step2：计算单个序列中相同单词片段的数量：

Seq1中相同word的位置：GC（1，5，7） CA（2，10） AT（3） TG（4） CG（6） CC（8，9） AC（11）

Seq2中相同word的位置：GA（1） AA（2，11） AT（3） TG（4） GC（5，7） CC（8） CA（10） AC（12）

Step3：随后将Seq1和Seq2中每一个单词出现的数量两两相乘，得到一个wordlist（Seq1存在而Seq2不存在的单词数量为0）：

Common words：GC=3*2=6；CA=2*1=2；AT=1*1=1；TG=1*1=1；CG=1*0=0；CC=2*1=2；AC=1*1=1

Step4：将所得结果求和，total=6+2+1+1+0+2+1=13

FASTA算法的数学公式为：Σw=|Lw（I）|*|Lw（J）|，其中I，J表示两个序列，Lw表示序列每个单词出现的次数，算法的复杂度为O（N1+N2），为了增加运行速度，可以将同一条对角线中的临近单词合并成一个。FASTA算法对DNA序列的搜索较于蛋白质更为灵敏，由于仅寻找最佳匹配，则有可能会遗漏掉其他有意义比对。

BLAST快速算法是另一种启发式局部算法[6][7]，目前被广泛应用并延伸出多种衍变算法，诸如：BLAST2，PSI-BLAST等等。基本思想是通过words的得分矩阵，取得质量更高的词表，然后在该词领域不断延伸寻找高分片段，当得分低于设定值时，中止寻找流程。

算法的主要步骤是：

Step1：定长分割序列，一般DNA序列的单词长度为3，蛋白质的单词长度为11

Step2：利用得分矩阵取出分值>值域的单词并组成单词表

Step3：两端延伸直至分值

Step4：将获得的高分片段相近合并

利用BLAST算法可以获得多条高分片段，能够保留更多有意义的匹配。

参考文献：

[1]Sniedovich， M. Dijkstra's algorithm revisited： the dynamic programming connexion[J]. Journal of Control and Cybernetics，2006，35（3）：599-620.

[2]Waterman， M. S. and M. Eggert . A new algorithm for best subsequence alignments with application to tRNA-rRNA comparisons[J]. J. Mol. Biol， 1987，197：723-728.

[3]Byung-Jun Yoon.Hidden Markov Models and their Applications in Biological Sequence Analysis[J]. PubMed Central， Curr Genomics，2009，10（6）：402-415.

[4]Krogh，A. An introduction to hidden Markov models for biological sequences[J]. Computational Methods in Molecular Biology . Elsevier， Amsterdam， The Netherlands， pp，1998：46-63.

[5]Lipman， DJ； Pearson， WR. Rapid and sensitive protein similarity searche[J]. Science， 1985，227（4693）：1435-41.

[6]Camacho， C.； Coulouris， G.； Avagyan， V.； Ma， N.； Papadopoulos， J.； Bealer， K.； Madden， T. L. BLAST+： Architecture and applications[J]. BMC Bioinformatics，2009，10：421.

[7]Kent， W. James. BLAT-The BLAST-Like Alignment Tool[J]. Genome Research，2002，12（4）：656-664.