阅读说明：本技术 一种用于基因比对的加速方法、装置、存储介质与服务器 (Acceleration method and device for gene comparison, storage medium and server ) 是由 谭光明 王元戎 张中海 于 2020-05-22 设计创作，主要内容包括：本发明提出一种用于基因比对的加速方法,其特征在于,包括：步骤1,获取待比对的基因序列,并将该基因序列中的各碱基存储于一压缩式BWT数据结构,得到一BWT序列,于该BWT序列中插入一计数点,得到一压缩序列；步骤2,获取一目标碱基,将该压缩序列中的该计数点以表项的形式存储于一查找表中,通过访存该查找表中的该计数点,对该压缩序列中的每一个碱基与所述目标碱基进行比对,得到比对结果；其中,所述压缩式BWT数据结构采用二进制位存储所述基因序列中的各碱基。(The invention provides an acceleration method for gene comparison, which is characterized by comprising the following steps: step 1, obtaining a gene sequence to be compared, storing each base in the gene sequence in a compressed BWT data structure to obtain a BWT sequence, and inserting a counting point in the BWT sequence to obtain a compressed sequence; step 2, obtaining a target base, storing the counting point in the compressed sequence in a lookup table in a form of table entry, and comparing each base in the compressed sequence with the target base through accessing the counting point in the lookup table to obtain a comparison result; wherein the compressed BWT data structure stores bases in the gene sequence in binary bits.)

一种用于基因比对的加速方法、装置、存储介质与服务器

技术领域

本发明涉及高性能计算领域，特别涉及一种面向基因比对领域的专用加速器及其加速方法。

背景技术

基因测序成本的下降加速了测序数据的产出，为数据处理带来了性能上的挑战。当前基因测序成本已降至1千美元以下，并且在未来还会继续大幅下降，这使得基因测序更加平民化，也导致测序数据***式增长，这为基因序列比对带来了极大压力。当下，基因序列比对的速度远远跟不上测序仪产生基因测序数据的速度，Illumina HiSeq2500测序仪产生长度为几百碱基的基因序列的吞吐率是每天150亿碱基(1500Gbps/day)，利用BWA-MEM比对工具的多线程版本，在双路24核心CPU服务器上试图将15GB大小的已测序片段比对到人类基因组参考序列之上，其测试结果最优性能为27Gbps/day，比测序仪吞吐速率慢50倍。因此，面对测序数据产出速度的不断提升，基因比对也必须不断提升其速率才能满足市场需求，采用定制结构对基因测序进行加速成为一个极受关注的方向。

基因比对应用以频繁的无规则片外访存为特征，在现有通用处理器结构上效率较低。如何用结构定制的方式设计一种面向基因比对领域的专用加速器，能够兼顾计算性能和能耗效率，成为一个难点。基因比对细粒度的应用特征，以及一系列与并发结构设计还不相适应的数据结构，都给充分发掘和利用应用的内在并行性带来了不小的难题。

发明内容

针对现有技术的不足，本发明提出一种面向基因比对领域的专用加速方法和加速装置，基于细粒度应用特征，提出了一系列与并发结构设计相适应的数据结构，通过数据通路和内存系统的定制来充分发掘和利用应用的内在并行性。

针对现有技术的不足，本发明提出一种用于基因比对的加速方法，其特征在于，包括：

步骤1，获取待比对的基因序列，并将该基因序列中的各碱基存储于一压缩式BWT数据结构，得到一BWT序列，于该BWT序列中***一计数点，得到一压缩序列；

步骤2，获取一目标碱基，将该压缩序列中的该计数点以表项的形式存储于一查找表中，通过访存该查找表中的该计数点，对该压缩序列中的每一个碱基与所述目标碱基进行比对，得到比对结果；

其中，所述压缩式BWT数据结构采用二进制位存储所述基因序列中的各碱基。

上述加速方法，其中，所述压缩式BWT数据结构中，所述碱基的存储方式为按位存储，分别用二进制编码表示A、C、G、T四种碱基。

上述加速方法，其中，其特征在于，所述A、C、G、T四种碱基分别采用二进制编码“00”、“01”、“10”、“11”表示。

上述加速方法，其中，所述步骤1进一步包括：于所述BWT序列中每隔N个碱基***一计数点，每N个该碱基与其前面邻近的该计数点组成一BWT区块，其中N＝16。

上述加速方法，其中，所述BWT区块中的该计数点包括4个64位整数，分别表示A、C、G和T四种碱基的频数。

上述加速方法，其中，所述查找表包括一索引，该索引的大小为8比特。

上述加速方法，其中，所述查找表包括256个表项，其中，每一所述表项的大小为32比特，所述表项用于记录该索引中A、C、G、T四种碱基的频数。

上述加速方法，其中，所述表项的的0-7比特表示碱基“A”的频数，8-15比特表示碱基“C”的频数，16-23比特表示碱基“G”的频数，24-31比特表示碱基“T”的频数。

上述加速方法，其中，所述步骤2进一步包括初始化所述查找表步骤：统计该查找表中每一所述表项的索引中包含的A、C、G和T四种碱基的数目，确定该表项的值。

上述加速方法，其中，所述压缩序列中每一个碱基与所述目标碱基进行比对的步骤进一步包括：

步骤21，将所述压缩序列中的每2个比特组成的比特对分别按位与所述目标碱基进行异或运算，得到第一比对结果；

步骤22，对所述第一比对结果进行或非运算，得到第二比对结果，该第二比对结果包括m个比特位，该m个比特位分别标识着所述压缩序列的m个碱基是否与所述目标碱基相同，其中m为正整数；

步骤23，计算所述第二比对结果中包含的“1”的数目M，M表示所述压缩序列中目标碱基的数目。

为实现本发明的另一目的，本发明还提供一种用于基因比对的加速装置，其特征在于，包括：

多个处理单元，用于对待比对的基因数据进行处理，且每一所述处理单元又包括一计算模块以及一片上存储模块；

一高带宽存储器，用于数据处理过程中的数据存取；以及

一二维网格片上网络，用于连接所述多个处理单元以及所述高带宽存储器；

其中，每一所述计算模块还包括：

一第一子模块，通过按位操作的LFMapBit映射算法实现，用于计算A、C、G、T四种碱基中的一种碱基的频数；以及

一第二子模块，通过基于快速查表的LFMapBitVec映射算法实现，用于计算A、C、G、T四种碱基的频数。

上述加速装置，其中，所述高带宽存储器为一多通道存储器(Multi-ChannelDRAM)。

上述加速装置，其中，所述片上存储模块包括：

一片上高速暂存存储器，用于存放所述查找表；

一FIFO存储器，用于缓存为所述处理单元预取的n条read数据；以及

一寄存器堆，该寄存器堆包括：

一常数向量寄存器，用于存储一计数数组；

一BWT区块偏移计算器，用于根据输入的BWT序列的所述索引，计算出所述索引所在的BWT区块在所述BWT序列中的偏移位置；

一计数点寄存器，用于存储所述BWT区块中的计数点；以及

一字向量寄存器，用于存储所述BWT序列中的字符串。

上述加速装置，其中，所述FIFO存储器缓存所述read数据的数目n＝2。

上述加速装置，其中，所述第一子模块包括一组合映射(Mask)阵列和一计数子模块，其中，所述组合映射(Mask)阵列又包括第一子阵列和第二子阵列。

上述加速装置，其中，所述第一子阵列用于获取一BWT区块和一目标碱基，并对该BWT区块中的碱基进行映射，同时并发地对该每一碱基与所述目标碱基进行异或运算，得到一比对向量；所述计数子模块用于获取所述比对向量，分别计算该比对向量中“0”和“1”的数量。

上述加速装置，其中，所述计数子模块为一种群计数(Popcount)归约树，该种群计数(Popcount)归约树为归约并发模式。

上述加速装置，其中，所述片上存储模块包括8个所述片上高速暂存存储器，该每一片上高速暂存存储器包括四个独立的读端口和四个存储体，所述片上高速暂存存储器将数据按顺序或交错地散列在该四个存储体上。

上述加速装置，其中，所述第二子模块通过所述片上高速暂存存储器的所述四个读端口查询同一所述查找表，并通过归约树计算所述BWT序列字符串中所有四种碱基的数目。

为实现本发明的另一目的，本发明还提供一种计算机存储介质，用于存储上述任一加速方法的计算机程序。

为实现本发明的另一目的，本发明还提供一种服务器，包括上述任一所述的加速装置。

由以上方案可知，本发明的优点在于：通过对基因比对核心算法进行并行化改进，使结构定制的优势能够更好发挥；结合应用的细粒度特征，对应用的核心操作模块进行数据通路和片上存储结构的定制，使得处理单元更加高效。

本发明设计一种面向基因比对领域的专用加速器，通过对基因比对核心算法进行并行化改进，使结构定制的优势能够更好发挥；结合应用的细粒度特征，对应用的核心操作模块进行数据通路和片上存储结构的定制，设计高效的处理单元。

以下结合附图和具体实施例对本发明进行详细描述，但不作为对本发明的限定。

附图说明

图1为本发明一实施例中的加速装置总体结构图；

图2为本发明一实施例中的按位表示的BWT数据结构图；

图3为本发明一实施例中的第一子模块微架构图；

图4为本发明一实施例中的第二子模块微架构图；

图5为本发明一实施例中的一个4端口的片上SPM的结构示意图；

图6为本发明一实施例中的局部碱基频数查找表；

图7(a)为本发明一实施例中的加速装置的在计数方面的总体性能展示图

图7(b)为本发明一实施例中的加速装置的在判定方面的总体性能展示图；

图8(a)为本发明一实施例中的加速装置在计数方面的片外带宽利用率比较图；

图8(b)为本发明一实施例中的加速装置在判定方面的片外带宽利用率比较图。

具体实施方式

为让本发明的上述特征和效果能阐述的更明确易懂，下文特举实施例，并配合说明书附图作详细说明如下。

如图1所示的本发明一实施例中的基因比对加速装置总体结构图，本发明中的基因比对加速装置主要包括加速装置专用处理单元(PE)、二维网格片上网络(On-chipNetwork)和多通道内存MCDRAM(Multi-Channel DRAM)三个部分。加速装置集成了众多PE单元并通过2D网格片上网络将他们互连起来，片上网络还连接了多通道内存MCDRAM(Multi-Channel DRAM)系统作为内存，这是一种新型3D结构的高带宽内存(HBM)。

加速装置中发处理单元(PE)的处理部件可以粗略分为计算部件和片上存储部件两类，片上存储部件主要包括寄存器堆、片上高速暂存存储器(Scratchpad Memory，SPM)以及缓存read数据的FIFO；处理单元(PE)的计算部件(Compute Units)主要包括两大模块：第一子模块和第二子模块，它们分别对应的是改进后的LFMapBit算法和LFMapBitVec算法的实现。其中，第二子模块，也叫LFMapBitVec模块可以利用查找表同时得到四种碱基的频数，用于实现前向计数过程第一子模块则通过专用位运算结构来计算某一种碱基的频数，可以用于后向计数和位置判定过程。

进一步地，在本发明的基因比对加速装置中，片上存储模块中包括一片上高速暂存存储器SPM，用于存放所述查找表；一FIFO存储器，用于缓存所述处理单元预取的n条read数据；以及一寄存器堆，该寄存器堆包括：一常数向量寄存器，用于存储一计数数组；一计数点寄存器，用于存储所述BWT区块中的计数点；一字向量寄存器，用于存储所述BWT序列中的字符串；以及至少8个通用寄存器，用于存储所述计算模块在计算过程中产生的其他变量和返回值。

更进一步地，为更好实现本发明的加速装置的良好性能，本发明还对加速装置中所涉及的软件模块进行了创新改进。如图2所示的本发明一实施例中的按位表示的BWT数据结构图，BWT(Burrows–Wheeler_transform)的较大体量加之其他辅助数据结构也会占用较大内存空间，为了减小内存的存储压力，提高访存性能，本实施例采用压缩式的BWT数据结构来对BWT序列进行存储。在这种压缩式数据结构中，A、C、G、T这四种碱基中的每一种碱基都由原来的字符存储(占一个字节)替代为按位存储，由于基因比对中只有这四种碱基，因此对于这四种碱基的编码只需要两位二进制数。本实施例把A、C、G、T四种碱基分别用“00”、“01”、“10”、“11”这四个两位二进制来编码，每个32位的整数代表着16个碱基。

根据FM-index算法的执行过程，本实施例在计算局部频数时需要对局部的BWT序列进行读取，而在这一步之前，还需要读取它前面邻近的计数点的数据(SegOA)，为了达到连续地址访存的效果，本实施例将计数点整合在BWT序列之中。在原BWT序列中每隔16个碱基***一个计数点，本实施例将每16个碱基以及其前面邻近的计数点组成的整段区域称为一个BWT区块(Rank)。其中Rank计数点由4个64位整数组成，分别表示A、C、C、T四种碱基的频数，大小为4*64bit，而Rank的序列区域大小为16*2bit＝32bit。在实际实现当中，本实施例每隔128个碱基***一个计数点，这样一个Rank的碱基序列区域大小为128*2bit＝32B，加之计数点的大小32B，整个Rank大小是32B+32B＝64B，正好是现有的基于Cache系统的处理器中一个Cache块的大小。数据结构存储优化：每个32位的整数代表着16个碱基，因为每个碱基是由两位二进制进行编码的。对局部的BWT序列进行读取时，读取它前面邻近的计数点的数据(SegOA)，为了达到连续地址访存的效果，将计数点整合在BWT序列之中。

加速装置的总体结构基于read处理的天然特征，每个处理单元(PE)可以独立的对read流进行批量处理，因此，各处理单元之间没有通信。设计所遵循的主要原则就是：第一，定制单个处理单元内“最短”数据通路，使单个PE能够达到极大吞吐率，从而使单个PE的带宽消耗能力提升；第二，通过PE数目扩展，充分支持数据并行，从而利用并行的内存访问进一步提升带宽利用效率。

对于本发明加速装置中的第一子模块，如图3所示为本发明一实施例中的第一子模块微架构图。第一子模块由组合映射Mask阵列和种群计数Popcount归约树两个主要部分组成，Mask阵列主要组成元素是第一子阵列和第二子阵列，其中第一子阵列可以为并发XOR阵列，第二子阵列可以为并发NOR阵列，其中XOR阵列属于“Map”模式，其输入为一个256位的BWT Rank字符串(看做一个“字”)和一个2位的目标碱基，需要进行“Map”操作的就是这个256位向量表示的128个碱基(每两位表示一个碱基)，需要同时并发的对每一个碱基与目标碱基进行对比(异或)。整体来说，Mask过程，即通过与目标碱基的位运算，将待搜索的“字”中的每个碱基转化为“0”或“1”信号，其中“0”表示该碱基不是目标碱基，“1”表示该碱基与目标碱基相同，需要计入。Popcount部分在Mask之后，实现“种群计数”功能，Popcount属于“Reduce”并发模式，128个位数值通过7级归约树求和，得到一个总和。其中，Popcount过程，即将转化后的“0/1”信号看做数字进行累加，从而计算出其中“1”的数目，也就原来的“字”中含有目标碱基的数目。其输入是Mask输出的128位向量，其中每一位是“0”或“1”，代表着原字符串中128个碱基是否与目标碱基相同。

对于本发明加速装置中的第二子模块，如图4所示为本发明一实施例中的第二子模块微架构图。第二子模块的作用是计算整个Rank字符串中所有四种碱基的数目，整个Rank的长度是256bit，存储在BWT字向量寄存器中。用8个高速暂存存储器SPM分别存储8张相同的查找表(Occurance)，每个高速暂存存储器SPM有4个独立的读端口，将Rank中的这32个字按顺序分为8组，每一组的4个字通过高速暂存存储器SPM的四个端口查询同一张表。在查表之后进行归约操作，用通用归约树实现，归约树(Reduce Tree)只需要支持加法操作即可。归约树的高度由归约的操作数的数目决定，对于32:1的归约，采用5层加法树。图5所示为一个4端口的片上高速暂存存储器SPM的结构示意图，由于对查找表的访问地址是无规则的，又要保证能同时支持四个独立的访问请求，我们将对高速暂存存储器SPM的设计采用内部分组的存储方式，即在一个高速暂存存储器SPM内部划分成四个存储体(Bank)，将数据按顺序或交错的模式散列在四个存储体上，以支持4个同时查表请求。

在本发明的另一实施例中，还提供了一种用于基因比对的加速方法，包括以下步骤：

步骤1、获取待比对的基因序列，采用一压缩式BWT数据结构存储该基因序列，得到一BWT序列，在该BWT序列中整合***一计数点，得到一压缩序列，其中，所述压缩式BWT数据结构采用二进制位存储所述基因序列中的碱基；

步骤2、将该压缩序列中的该计数点以表项的形式存储在一查找表中，获取目标碱基，通过访存该查找表中的该计数点，对该压缩序列中每一个碱基与目标碱基进行比对，得到比对结果。

在上述加速方法中涉及LF映射算法的改进：最原始的LF映射的实现方式，其中的碱基是按字节存储，并且再进行碱基局部频数的计算时，对碱基序列的扫描也是以一个字节为单位进行的。为了提升存储效率，将数据结构的存储方式改进为按二进制位来存储碱基。同时，采用快速查表的方法用访存替代计算，将所有四种碱基的局部频数都以表项的形式存储在“查找表”中。

按位操作的LF映射算法LFMapBit，具体来讲，要计算一个“字”(假设一个字为一个字节，占用8比特)中某个碱基的数目，需要以下两个步骤：第一，将这个字(8比特)中的每2比特(也就是一个碱基)组成的“比特对”分别按位与目标碱基(2比特)进行异或运算，这时再将每个“比特对”得到结果(也是“比特对”)进行“或非”运算，最终，原来的一个“字”会变为4个比特位，分别标识着这个字中的四个碱基是否与目标碱基相同(1表示相同，0表示不同)。第二，字中目标碱基的数目就是上一步结果中“1”的数目，因此，我们只需要统计上面结果中“1”的个数，也就是简单地将得到的4个比特标识看做数字进行加法运算即可。

基于快速查表的LF映射算法LFMapBitVec，具体来讲，我们将8比特大小的一个“字”作为查找表的索引，这也决定了查找表的表项有256个，如图6所示的局部碱基频数查找表，每个表项大小是32比特(4字节)，记录着这个8比特“字”中A、C、G、T四种碱基的频数，其中0-7比特表示碱基“A”的频数，8-15比特表示碱基“C”的频数，16-23比特表示碱基“G”的频数，24-31比特表示碱基“T”的频数。表的初始化是通过统计每个表项的索引中包含的四种碱基的数目来确定表项的值，例如，对索引为88的表项初始化，由于“88”的二进制表示为“01011000”，将它与A(00)、C(01)、G(10)、T(11)四种碱基进行对比，发现其中包含1个“A”、2个“C”、1个“G”和0个“T”，因此(1,2,1,0)便是此表项的值。

在本发明的一实施例中，本发明的加速装置的数据处理流程具体为：对于输入的fastq文件，加速器的每个处理单元(PE)可以独立的对read流进行批量处理。对于流式的read数据，其存储地址是连续的，加速装置采用预取策略对后续的read进行预取，并缓存在FIFO中。假设缓存read的FIFO可以容纳的read数目为N，实验发现，N＝2对于本发明的加速装置已经足够，因为对于当前常用的100bp左右长度的read来说，计算部件处理完一条read通常需要上百次LF映射迭代，因此预取数目没有必要很多。计数数组AC存放于256位的常数向量寄存器中。

在每次LF迭代过程中，目标碱基选定后，Rank偏移计算器会根据输入的BWT索引(也就是sp或者ep)计算出索引所在的Rank在BWT序列中偏移位置，计算出地址后，相应的Rank数据就会整个从内存中取出来，Rank计数点存放在256位的计数点寄存器(CheckpointRegister)中，同时BWT字符串存放在256位的字向量寄存器中(BWT Word Register)，此外，本发明的加速装置有8个64位通用寄存器(GPRs)来暂存计算过程中的其他变量和返回值等。片上高速暂存存储器是用来存放频数查找表的。

在本发明的另一实施例中，还提供了一种计算机存储介质，用于存储上述任一加速方法的计算机程序。

在本发明的另一实施例中，还提供了一种服务器，包括上述任一所述的加速装置。

综合所述，本发明主要包括以下方面：

1.算法-结构协同设计。从算法与结构协同设计的角度对应用特征进行了深度的剖析，找到与硬件结构设计的契合点，通过关键算法的并行化改进。

2.基于对比对应用在传统处理器结构上的性能瓶颈的观察和分析，提出了一种基因序列比对加速器结构，其主要设计原则是充分利用应用特有的细粒度并行特征并充分挖掘访存并行性，设计定制的专用数据通路和大规模数据并行结构，从而高效利用现有内存提供的访存带宽。

为了更好体现本发明加速装置性能，本发明还为其配备了两种软件实现的对比基准(Baseline)，其运行平台和软件配置不同，第一个基准是在E5-2670CPU平台上的单线程软件实现，其处理器核心运行在2.6GHz的主频上，其内存系统为32GB容量DDR3-1333内存，利用4个内存通道，可以提供40+GB/s的内存带宽；第二个基准为Intel KNL平台上的256线程的软件实现，KNL众核协处理器平台有68个乱序执行核心，其处理器核心运行主频为1.6GHz，KNL使用基于多通道动态随机访存MCDRAM的高带宽内存(HBM)，这种内存可以提供高达400GB/s的外部内存带宽。对于Read数据集，本发明采用来自Illumina platinumgenomes的101bp长度的真实人类基因片段和用Wgsim片段产生器构造的大于200bp长度的模拟Read数据。

图7(a)和图7(b)展示了本发明的加速装置相对于两种通用处理器平台上的最优软件实现的归一化加速比，其中Xeon E5基准平台的性能被归一化为“1”。本发明评测了以人类全基因组以及各染色体为参考序列时的加速比，其加速比差异主要来源于各场景下参考序列的大小差异，参考序列体量大，其在通用平台上cache失效率越高，如全基因组参考序列，因此相比之下其加速比更高。

结果显示，本发明的加速装置相对于Xeon E5 CPU处理器上的单线程最优软件实现达到了720倍至990倍的加速比，而相对于KNL处理器上的256线程配置的软件实现也获得了19倍至44倍的加速比。本发明评测了以人类全基因组以及各染色体为参考序列时的加速比，其加速比差异主要来源于各场景下参考序列的大小差异，参考序列体量大，其在通用平台上cache失效率越高，如全基因组参考序列，因此相比之下其加速比更高。结果显示，本发明的加速装置相对于XeonE5 CPU处理器上的单线程最优软件实现达到了720倍至990倍的加速比，而相对于KNL处理器上的256线程配置的软件实现也获得了19倍至44倍的加速比。

图8(a)和图8(b)展示了本发明的加速装置相对于KNL平台上的最优软件实现的归一化片外内存带宽利用率，其中KNL的带宽利用率被归一化为“1”。结果显示，相比于KNL平台上256线程配置的最优软件实现，本发明的加速装置的内存带宽利用率是前者的24倍至86倍。

当然，本发明还可有其它多种实施例，在不背离本发明精神及其实质的情况下，熟悉本领域的技术人员当可根据本发明作出各种相应的改变和变形，但这些相应的改变和变形都应属于本发明所附的权利要求的保护范围。