基于HLS的BWT算法的硬件加速方法

2021-03-11 03:34:54馬毅超張曉盼王鑫

電子技術與軟件工程 2021年21期

關鍵詞：優化

馬毅超張曉盼王鑫

（陜西科技大學電子信息與人工智能學院陜西省西安市 710021）

1 引言

序列對比[1]是指將測序儀單次測序所得到的長度較短的堿基序列（reads）映射到參考基因組上的過程。序列對比經常被用于分析序列片段的相似性，通過對比結果得出物種之間的同源性及其生物學特征，以此來判斷生物的健康狀況以及進化信息等。近年來，生物學數據庫中的序列數據呈現爆炸式增長，以美國國家生物技術信息中心(NationalCenterforBiotechnologyInformation，NCBI)建立的DNA 序列數據庫GenBank 為例，DNA、RNA 和蛋白質的序列數據大約每6 個月就會翻一番[2]。新一代測序技術也將產生海量的測序數據序列（reads），生物信息學將面臨著對這些海量數據處理的嚴峻挑戰。因此，通過研究序列對比問題可以加快對數據的處理速度，從而應對生物信息學對大量數據研究的需求。

M.Burrows and D.J.Wheeler 于1994年提出了BWT 算法[3]，該思想首先被用于數據的壓縮算法[4]，之后逐漸被應用于字符串的匹配。BWT 算法于2007年首次應用于DNA 序列對比問題中，2009年以后普遍用于將DNA 短序列映射到參考基因組。如Bowtie，BWA，SOAP2等主流的序列對比軟件都采用了BWT技術[2]，但是BWT 算法運算量大，當數據規模較大時，從參考基因組中查找短序列reads 的時間較長。因此，本文提出了一種基于HLS 的并行方法對BWT 算法進行加速。

2 Vivado High-level Synthesis介紹

Xilinx 公司推出的Vivado High-Level Synthesis (HLS)工具[5]可以將C 或者C++規范轉換為寄存器傳輸級（RTL）實現，并將其合成為現場可編程門陣列（FPGA）。FPGA 則提供了一中大規模的并行架構，它在性能、成本和功率方面優于傳統處理器。近年來，高級綜合（HLS）工具不斷得到改進，同時提供了各種硬件優化的可能性。使用不同的優化及其組合會導致同一內核的不同版本，就像使用不同的編譯器標志時從同一源代碼生成不同的二進制文件一樣，這些實現可能在性能/面積/功率方面有所不同。由于其相較于傳統的FPGA 開發周期短，因此越來越多的被應用于各類算法的加速優化[6-11]。其開發流程如圖1所示。High-Level Synthesis 包含以下階段：

圖1

（1）時序安排。根據以下內容確定每個時鐘周期發生的操作：

①時鐘周期或者時鐘頻率的長度；

②在目標芯片上完成該操作所需要的時間；

③用戶指定的優化指令；

（2）綁定。確定實現每個計劃操作所需要的硬件資源。

（3）控制邏輯提取。提取控制邏輯以創建有限狀態機（FSM），該狀態機對RTL 中設計中的操作進行排序。

3 算法原理

3.1 后綴數組[12]

3.2 Burrows-Wheeler Transform[13]

圖2

圖3

3.3 Backward Search[13]

在下述例子中應用了BWT 的first-last 屬性，介紹使用BWT的后向搜索在字符串S=CTATATGC 精確匹配其子串P=TAT，匹配過程如圖4所示，算法的偽代碼如圖5所示，得到reads 所在的區間[sp, ep]后再由SA[i]可以還原reads 在原始字符串S 中的位置，其中根據圖4結果有，在圖3中有SA[6]=1，SA[7]=3，找到子串P 在S 中位置為1 和3。

圖4：Backward Search

圖5：Backward Search 偽代碼

從偽代碼中可以看出函數C(symbol) 和Occsymbol(L,i) 在Backward Search 精確匹配算法和LFsymbol(i)算法中都被反復調用。同時由于函數Occsymbol(L,i)計算量大，耗費了大量的時間在該計算過程中。因此對這部分代碼進行優化，可以提高運算速度，從而可以有效提高DNA 序列比對的效率。

4 BWT的HLS實現

實驗中使用的參考基因組是從NCBI 官網（https://www.ncbi.nlm.nih.gov/nuccore/NC_008253.1?report=fasta）下載的Escherichia coli 536, complete sequence，序列編號為NC_008253.1，數據大小約為5MB。短序列reads 總大小為1MB,每條reads 長70bp,共5000 條。

4.1 串行方案

針對在算法原理介紹部分提出的Occsymbol(L,i)函數計算量大并且被頻繁調用的問題,在代碼設計時,提前將計算結果用數組保存起來,以便后續查找。HLS 實現BWT 算法加速流程如圖6所示，其中每一條reads 都需要通過Test Bench 傳入Search 函數來得到匹配結果，函數Search 自然成為高級綜合階段的頂層函數。在Search中首先需要計算在后續匹配過程中使用到的BWT, SA ,Occ 數組，這三個數組只和5M 的參考基因組相關，在之后所有的匹配中只與這三個數組相關，與原來的參考基因組無關，因此只需在第一條reads 匹配時創建，同時也就意味著第一條reads 的匹配時間要遠大于其他各條reads 的匹配時間。

表1：串行并行方案耗時比較

表2：串行并行方案所耗資源比較

圖6：硬件加速流程

圖7

4.2 并行方案

4.2.1 UNROLL 優化

針對創建BWT，SA，Occ 三個數組時耗費較長時間這一現象，通過HLS 對代碼中的大量循環（部分代碼如7 所示）添加unroll等指令對程序進行優化。對于圖8(a)中的for 循環生成的RTL 邏輯電路執行過程如圖8(b)所示，每次循環使用的時同一套加法器，每次循環分時復用該加法器。

圖8

對for 循環添加unroll 指令如圖8(a)，該指令通過設置參數factor=3 將上述的循環電路復制了3 份，這些復制出的電路可以同時處理不同數據，處理過程如圖8(b)。由于該優化指令要求循環的邊界是常量，因此當循環邊界是變量時，需要在循環前設置斷言，來設置循環的最大次數。

4.2.1 ALLOCATION 優化

在Backward Search 中每次匹配一個字符都需要使用LF 函數不斷更新reads 可能存在的位置區間[sp，ep]，在計算sp 和ep 時程序串行執行，調用兩次LF 函數，執行過程如圖9(a)，通過HLS 對程序添加ALLOCATION 指令如圖9(b)，將LF 函數實例化為兩個module 如圖9(c)，可并行得出sp 和ep 的值。

圖9

4.3 資源優化

在提升算法效率的同時也要盡可能的減少FPGA 資源的使用。若將全量的索引數據全部存入RAM 將會對系統造成相當大的負擔。為了減少資源占用，計數數組num_symbol 均以32 為閾值進行壓縮。

Occsymbol(L,i)表示最后一列L 中的前i 個位置中字符symbol 出現的次數,由于每次計算Occ 都要遍歷數組L 將耗費大量的時間,因此使用計數數組num_symbol 加快Occ 的計算。圖11 中，每個計數數組num_symbol 及其對應的L 數組的位置范圍被標記在同一陰影級別。代碼實現如圖12 所示,其中num_symbol[i]數組中存儲的數組L 中下標0~(i+1)*32-1 之間字符symbol 的個數，右側在計算Occ 時將不能被32 整除的余數部分在數組L 中進行遍歷操作，遍歷時間也會隨著閾值的增大而增大，因此閾值不宜過大，最后再加上可以被32 整除部分的計數值就可以得出結果。將在壓縮計數數組和計算Occ 時使用到的乘法運算符及除法運算符全部替換為“＜＜”(左移運算符)和“＞＞”（右移運算符），減少DSP 資源的消耗，同時可以提高運算效率。

圖11：L 與num_symbol 的映射關系

圖12：計數數組生成及Occ 計算

4.4 仿真結果

將仿真波形圖在Vivado 中打開，如圖13 所示，其中ap_clk 為時鐘信號，ap_rst 為復位信號。n 為每條reads 長度長度為71（包含字符串結束符’