并行计算在生物序列比对中的应用

郝静 刘雅坤

摘 要：生物信息学作为一门综合运用分子生物学、数学和计算机等学科的理论和方法的交叉学科为阐明和理解海量DNA序列数据所包含的生物意义提供了可能。而海量的数据信息，使得生物序列比对计算需要耗费大量的时间，本文中针对此问题，实现Windows API，OpenMP的并行算法，比较各自的加速比和最优比对序列。

关键词：生物序列比对 Windows API OpenMP

一、问题背景与提出

生物信息学最首要的任务就是从大量的生物信息数据中，提取有生物价值的信息，而序列比对是当前最重要，最基本的方法，是指两个或多个序列按字母比较， 尽可能确切地反映它们之间的相似和相异性，阐明序列之间的同源关系。在序列分析中， 将未知序列同已知序列进行相似性比较是一种强有力的研究手段， 从序列的片段测定， 拼接， 基因的表达分析， 到 RNA 和蛋白质的结构功能预测， 物种亲缘树的构建都需要进行生物分子序列的相似性比较。

一般来说， 评价生物序列比对算法的标准有两个： 一为算法的运算速度， 二为获得最佳比对结果的敏感性或准确性。然而，随着测序技术的发展，生物序列分析要求新的软件方法和计算平台，在分析大规模序列数据和解决复杂问题上，并行计算发挥着巨大的作用。

二、模型建立

序列比对是运用某种特定的数学模型或算法， 找出两个或多个序列之间的最大匹配碱基或残基数， 比对算法的结果在很大程度上反映了序列之间的相似性程度以及它们的生物学特征。序列比对根据同时进行比对的序列数目多少可分为双序列比对和多序列比对，从比对范围考虑也可分为全局比对和局部比对。

2.1 两两序列比对

两两序列比对，就是把两条未知的序列进行排列，通过字母的匹配，删除和插入操作，使得两条序列达到同样长度，在操作的过程中，尽可能保持相同的字母对应在同一个位置。

通常用得分矩阵描述序列两两比对， 两条序列分别作为矩阵的两维， 矩阵点记录两个维上对应的两个DNA序列的相似性分数， 分数越高则说明两个DNA序列越相似。

2.2 序列比对的基本定义

在生物分子信息处理过程中， 将生物分子序列抽象为字符串， 其中的字符取自特定的字母表。例如，DNA 序列由四种核苷酸组成， 用“A”， “T”， “C”， “G”代表四种碱基， 其复杂度为 4。

生物序列比对可以看作字符串的比对，用如下的定义来描述序列的比对与相似性。

空位的引入是为了补偿由于变异而产生的插入或缺失，每引入一个空位，比对的分值都会有所扣除：

其中， 表示空位的总罚分； 表示初始化一个空位的罚分； 表示空位延伸一个间隔的罚分。

Def4 全局最优比对和局部最优比对

对于两条序列S和T的全局比对可以用序列 和 来表示，其中，

l） 和 分别是在S和T中加入一些空位而得到的序列；

2） ， 为序列 的长度。

将序列 和 相应的位置进行一一比较，其分值Score可用如下的公式来表示：

序列S和T的全局最优比对：指在S和T的所有比对中相似性分值Socer最大时所对应的比对。

（2）序列S和T的局部最优比对：用 来表示，其中 分别为S和T的子序列，且 的Score分值是S和T中所有子序列比对分值的最大值。

三、局部最优序列比对串行算法描述

3.1 基本思想

本文针对两两序列最优局部比对的Smith-Waterman方法，利用初始条件和迭代关系得到一个得分矩阵，选取其中得分最高的子序列的末端，然后利用动态规划的方法回溯寻找，得到局部最优比对序列。

对于两个长度分别为n和m 的DNA序列：

构造得分矩阵 ，用来存放所有可能的对比结果。

初始条件：D（i，0）=D（0，i）=0

递归关系：

其中， ， 分别为在序列S和T中添加一个长度为x，y 的空位的罚分。

在得分矩阵D中找出 ： ，则 表示序列S和T的局部最优比对分值。从 开始，按照从右下到

左上的方向，对矩阵D进行回溯就可以求得局部最优比对序列。

3.2 算法描述

为简化问题，本文中假设序列S和T中分别只可加入一个空位，且空位间隔为1，则设 。

从而递归关系简化如下：

四、某问题的并行算法描述

4.1 基于Windows API的多核并行算法的设计

Windows API多核并行算法利用WINAPI定义线程函数，然后主函数里创建线程，将线程函数导入创建的线程中运行。具体步骤如下：

（1）主函数中利用CreateThread函数创建线程；

（2）定义线程函数，两个线程分别调用，并根据线程号均分任务；

（3）为避免数据竞争，利用临界区比较中间得分矩阵的最大值，得到最终结果；

（4）利用回溯在（3）得到的得分矩阵中找到局部最优序列比对。

4.2 基于Open MP的多核并行算法的设计

OpenMP是一种面向共享内存以及分布式共享内存的多处理器多线程并行编程语言，其执行模式采用Fork-Join的形式，当程序执行时，主线程生成（Fork）一组线程，随着程序的执行，在并行结束后，这组线程汇合（Join）成主线程。具体步骤如下：

（1）利用threadprivate指令将全局变量thread_xl复制到各个线程中；

（2）使用#pragma omp parallel for 将已知序列S均分到两个线程中，实现循环并行化，分别得到两个得分矩阵；

（3）利用#pragma omp parallel实现并行，求各线程得分矩阵中的最大元素及其位置；

（4）利用原子操作对两个得分矩阵中的最大值进行比较，选其较大者；

利用回溯在（4）得到的得分矩阵中找到局部最优序列比对。

五、算法实现

5.1 Windows API

1.任务分配：

HANDLE hThread[numthread]；

int myid[numthread]；

InitializeCriticalSection（&cs）；

for（int i=0；i

{myid[i]=i；hThread[i]=CreateThread（NULL，0，API_getD，&myid[i]，0，NULL）；}

WaitForMultipleObjects（numthread，hThread，TRUE，INFINITE）；

DeleteCriticalSection（&cs）；

2.线程函数：

DWORD WINAPI API_getD（LPVOID arg） {

int threadID=*（（int*）arg）；

XL thread_xl={{0}，{0}，{0}，{0}}；

int p=1，n=5，m=6；

if（threadID==1）

{p=5；n=9；}

getD（p，n，m，thread_xl）；

D_Max（n，m，thread_xl）；

EnterCriticalSection（&cs）；

if（thread_xl.max>xl_2.max）

xl_2=thread_xl；

LeaveCriticalSection（&cs）；

return 0；}

5.2 OpenMP

1.并行计算得分矩阵

#pragma omp parallel for

for（int i=1；i<9；i++）{

for（int j=1；j<6；j++）{

thread_xl.score[i][j]=Score（s[i]，t[j]）；

int a=thread_xl.D[i-1][j-1]+thread_xl.score[i][j]；

int b=thread_xl.D[i-1][j]-W；

int c=thread_xl.D[i][j-1]-W；

thread_xl.D[i][j]=Max（a，b，c）；}}

2.计算最终得分矩阵中最大值以及位置

#pragma omp parallel

int myid=omp_get_thread_num（）；

D_Max（9，6，thread_xl）；

#pragma omp critical

if（thread_xl.max>xl_1.max）{xl_1=thread_xl；}

数值实验和结论

为了比较串行与多种并行算法的效率，以串长分别为9和6的序列为例，进行数值试验，得到结果如下。

以串行计算时间为基础，采用API，OpenMP算法使用两个核计算的加速比分别为：0.430，0.667；加速比不高，这是因为针对本例中的S序列和T序列，序列长度短，在寻找最优比对的过程中，计算量不大，所以串行也可以很快的完成计算；而并行在创建线程时消耗了大量时间，因此出现加速比小于1的情况。

得到的最优序列比对分别为：串行与OpenMP：tcggc，t-ggc；API：ggc。结果的不同是由于这几种并行计算算法的差异导致的。在具体分配任务时，API与MPI真正的将S序列分成两段，得到两个得分矩阵选择含有最大得分的矩阵进行回溯；而在OpenMP中，虽然将S序列分成两段，但是得到的两个得分矩阵是存放在一个大矩阵中进行回溯。

从以上结果的分析中，可以看出，OpenMP算法的并行设计更适合此类问题。

为了比较不同迭代次数下多种计算方式的计算效率，在程序实现中，增加随机赋值已知序列与未知序列部分，增加其序列长度，得到计算结果如上所示。总体上，三种并行方式的加速比随着迭代次数的增加是提高的，但受到计算量的限制，加速比仍未体现出并行计算在大计算量运算中的优势，但是较上例已经有所提高。

参考文献

[1] 陈华.高性能并行计算.中国石油大学（华东）[M].2016

[2] 张法.基于Smith_Waterman算法的并行分而治之生物序列比对算法[J].中国科学：技术科学.2004.34（2）：190-199

作者简介：

郝静（1994—），女，汉族，山东烟台人，中国石油大学（华东）理学院，2013级本科生，数学与应用数学专业。

刘雅坤（1995—），女，汉族，山东青岛人，中国石油大学（华东）理学院，2013级本科生，数学与应用数学专业。