李金光 唐友 李丹

基金項(xiàng)目：吉林省科技發(fā)展計(jì)劃項(xiàng)目；項(xiàng)目名稱(chēng)：基于數(shù)據(jù)挖掘技術(shù)的全基因組選擇方法研發(fā)及云計(jì)算平臺(tái)體系構(gòu)建；項(xiàng)目編號(hào)：YDZJ202201ZYTS692。

作者簡(jiǎn)介：李金光（1999— ），男，河北廊坊人，碩士研究生；研究方向：先進(jìn)控制技術(shù)與控制系統(tǒng)集成。

*通信作者：李丹（1980— ），男，吉林吉林人，副教授，博士研究生；研究方向：電氣工程。

摘要：近年來(lái)，多核技術(shù)在基因測(cè)序方面有著廣泛的應(yīng)用。因此，基因測(cè)序的并行化成為目前國(guó)內(nèi)熱點(diǎn)。在基因測(cè)序領(lǐng)域，隨著高通量測(cè)序技術(shù)的迅猛發(fā)展和應(yīng)用越加廣泛，其無(wú)時(shí)無(wú)刻不在產(chǎn)生大量的數(shù)據(jù)。對(duì)此，串行程序進(jìn)行單個(gè)程序運(yùn)行已經(jīng)無(wú)法滿(mǎn)足社會(huì)對(duì)基因測(cè)序的需求。為了同時(shí)分析數(shù)千個(gè)樣本，測(cè)序需要快速、靈活和內(nèi)存高效地實(shí)現(xiàn)。文章認(rèn)為，對(duì)一個(gè)名為ANGSD的程序進(jìn)行OpenMP的程序并行的軟件優(yōu)化設(shè)計(jì)和可行性研究非常必要。

關(guān)鍵詞：OpenMP；下一代測(cè)序；多線(xiàn)程

中圖分類(lèi)號(hào)：TP39 文獻(xiàn)標(biāo)志碼：A

0 引言

測(cè)序技術(shù)經(jīng)過(guò)多年發(fā)展擁有了3個(gè)階段成果，分別是第一代測(cè)序技術(shù)、下一代測(cè)序技術(shù)（NGS）和最新的第三代單分子測(cè)序技術(shù)。下一代測(cè)序技術(shù)（Next Generation Sequencing，NGS）又被叫作高通量測(cè)序技術(shù)，此技術(shù)一次測(cè)序可以檢測(cè)幾十萬(wàn)到幾百萬(wàn)條序列。由于此技術(shù)具有高通量、強(qiáng)穩(wěn)定、高準(zhǔn)確度等優(yōu)點(diǎn)，因此在商業(yè)上被廣泛應(yīng)用［1］。

現(xiàn)有的通用多樣本NGS分析程序的例子有單線(xiàn)程SAMtools （C）和多線(xiàn)程GATK（Java）以及ANGSD。這3個(gè)程序之間有許多不同之處，但ANGSD的關(guān)鍵優(yōu)勢(shì)在于：（1）允許與原始測(cè)序數(shù)據(jù)直接相關(guān)的多種輸入數(shù)據(jù)類(lèi)型（文本堆積、二進(jìn)制基因型可能性文件、VCF文件）；（2）允許用戶(hù)選擇多種方法中間分析（例如計(jì)算GL的不同方法）；（3）包括一組在任何其他軟件中都沒(méi)有實(shí)現(xiàn)的下游分析的實(shí)現(xiàn)。ANGSD中很多方法都是基于GLs。ANGSD支持4種不同的模型來(lái)計(jì)算GLs：（1）重新校準(zhǔn)的SOAPsnp模型；（2）原GATK模型；（3）SAMtools 1.16 +修改Maq模型；（4）特定類(lèi)型的錯(cuò)誤模型。這些GL模型中的排序錯(cuò)誤率要么是固定的（由qscores獲得），要么是從數(shù)據(jù)中估計(jì)出來(lái)的。4個(gè)實(shí)現(xiàn)的GL模型假設(shè)為二倍體樣本。工作流程分為2個(gè)步驟：（1）ANGSD生成特定的分析輸入數(shù)據(jù)。（2）根據(jù)ANGSD輸出，使用二級(jí)關(guān)聯(lián)程序進(jìn)行下游分析。對(duì)于簡(jiǎn)單的測(cè)試，例如ABBA-BABA/D-statistic，次要程序可以是簡(jiǎn)單的Rscript；對(duì)于計(jì)算密集型的方法，次要程序可以是多線(xiàn)程的c/c++程序［2］。OpenMP是多核系統(tǒng)最流行的并行編程之一。因此，本文就OpenMP并行在ANGSD軟件的應(yīng)用可行性做展開(kāi)論述。

1 研究意義

20世紀(jì)70年代，Walter Gilbert和Frederick Sanger發(fā)明了人類(lèi)史上第一臺(tái)測(cè)序儀，并通過(guò)測(cè)序儀探測(cè)到第一個(gè)基因組序列噬菌體X174，全長(zhǎng)擁有5 375個(gè)堿基。

在Sanger基因測(cè)序技術(shù)的催動(dòng)下，人們對(duì)生命的本質(zhì)產(chǎn)生極大的興趣。在此動(dòng)力下，人們對(duì)生命科學(xué)的研究逐漸步入了基因組學(xué)的時(shí)代。40多年來(lái)，人們通過(guò)不斷努力已將測(cè)序技術(shù)發(fā)展到了可觀(guān)的境界，已經(jīng)由第一代發(fā)展到了第三代測(cè)序技術(shù)。

由Sanger發(fā)明的最初的測(cè)序方法被人們稱(chēng)為第一代測(cè)序技術(shù)，雖然是第一代技術(shù)，但至今仍在被廣泛應(yīng)用。同時(shí)，它也有著無(wú)法掩蓋的缺陷。它進(jìn)行一次測(cè)序只能得到一條長(zhǎng)度為700～1 000個(gè)堿基的序列，這個(gè)缺陷成為在世界舞臺(tái)上更進(jìn)一步的掣肘，無(wú)法滿(mǎn)足現(xiàn)代科學(xué)發(fā)展所需要的龐大生物基因序列。

高通量測(cè)序（High-Throughput Sequencing，HTS）彌補(bǔ)了第一代基因測(cè)序的缺陷，運(yùn)行一次便可以同時(shí)獲得幾十萬(wàn)至幾百萬(wàn)條核酸分子的序列，因?yàn)樵赟anger基因測(cè)序技術(shù)的基礎(chǔ)上進(jìn)行了突破，所以高通量測(cè)序也被稱(chēng)為新一代測(cè)序 （Next Generation Sequencing，NGS）或第二代測(cè)序。

第二代測(cè)序技術(shù)雖然對(duì)各個(gè)領(lǐng)域產(chǎn)生了極大的影響，但是也有缺點(diǎn)。第二代測(cè)序技術(shù)獲得單條序列長(zhǎng)度很短，如果想得到更為準(zhǔn)確的基因序列信息，就不得不依賴(lài)于較高的測(cè)序覆蓋度和更準(zhǔn)確的序列拼接技術(shù)。由于技術(shù)問(wèn)題，最終得到的結(jié)果中會(huì)存在一定量的錯(cuò)誤信息，誤導(dǎo)研究人員的判斷。

為了防止結(jié)果中出現(xiàn)錯(cuò)誤信息擾亂人的認(rèn)知，科研人員在保證高通量測(cè)序的基礎(chǔ)上，發(fā)明了第三代測(cè)序技術(shù)，也稱(chēng)為單分子測(cè)序技術(shù)（Single-molecule sequencing technology）。該技術(shù)在上一代測(cè)序技術(shù)上做出了突破，解決了單條序列長(zhǎng)度過(guò)短的問(wèn)題，如今能夠直接得到長(zhǎng)度在數(shù)萬(wàn)個(gè)堿基的核酸序列信息。

目前，基因測(cè)序技術(shù)在包括胚胎的植入前遺傳學(xué)診斷研究、法醫(yī)學(xué)的少量DNA測(cè)序、表觀(guān)遺傳學(xué)和考古學(xué)、物種進(jìn)化演替過(guò)程等眾多領(lǐng)域得到廣泛應(yīng)用。

就當(dāng)前的市場(chǎng)勢(shì)態(tài)而言，第二代短讀長(zhǎng)序列檢測(cè)技術(shù)在世界測(cè)序技術(shù)市場(chǎng)上一直占據(jù)著絕對(duì)的優(yōu)勢(shì)地位，而第三代測(cè)序技術(shù)也在近年來(lái)實(shí)驗(yàn)中迅速發(fā)展。

若要進(jìn)行ANGSD軟件的簡(jiǎn)化提速，就需要對(duì)ANGSD軟件運(yùn)行的流程進(jìn)行分析：可以通過(guò)程序運(yùn)行來(lái)確定ANGSD讀取數(shù)據(jù)、運(yùn)算以及輸出數(shù)據(jù)的時(shí)間占比；還可以通過(guò)Pthread和OpenMP的兩種并行語(yǔ)言的對(duì)比進(jìn)行分析，以此來(lái)判斷ANGSD軟件是否擁有簡(jiǎn)化的可能性。

2 ANGSD數(shù)據(jù)處理流程

ANGSD是一個(gè)多線(xiàn)程程序。該程序可以計(jì)算各種匯總統(tǒng)計(jì)，利用下一代測(cè)序數(shù)據(jù)中的完整信息，通過(guò)直接處理原始測(cè)序數(shù)據(jù)或使用基因型可能性，進(jìn)行關(guān)聯(lián)映射和群體遺傳分析。該程序不僅支持多種輸入格式，包括BAM和beagle基因型概率文件，而且允許用戶(hù)在現(xiàn)有方法的組合之間進(jìn)行選擇，并可以執(zhí)行其他地方?jīng)]有實(shí)現(xiàn)的分析。此篇文章便對(duì)輸入為BAM基因型概率的文件進(jìn)行研究。

本文通過(guò)利用gdb調(diào)試工具，來(lái)判斷程序函數(shù)的調(diào)用關(guān)系，得出圖1為ANGSD軟件的函數(shù)調(diào)用流程。命令如下：

gdb ./angsd

set args-GL 1-doGlf 2-b bam2.filelist-doMajorMinor 1 -doMaf 1

b angsd.cpp：main

b shared.cpp：main_analysis

set print pretty on

r

以上已經(jīng)得出ANGSD軟件運(yùn)行處理的流程，接下來(lái)為了更直觀(guān)地表現(xiàn)出哪一部分在軟件中可以進(jìn)行簡(jiǎn)化更改，本文對(duì)程序選取不同規(guī)格的Bam文件進(jìn)行測(cè)序，觀(guān)測(cè)在虛擬機(jī)以及服務(wù)器的運(yùn)行速度。

從表1可以發(fā)現(xiàn)不同容量Bam文件分別在reads運(yùn)行時(shí)間占CPU總時(shí)間的5%，5%，7.5%，7.8%，所以可以得出結(jié)論讀取Bam文件reads記錄時(shí)間約為總時(shí)間的7%。

從表2可以發(fā)現(xiàn)不同容量Bam文件分別在輸出beagle基因型概率文件時(shí)間占CPU總時(shí)間的22%，20%，11.4%，10.9%，所以可以得出結(jié)論輸出beagle基因型概率文件約為總時(shí)間的10%。

從表3可以發(fā)現(xiàn)不同容量Bam文件分別在計(jì)算基因型概率時(shí)間占CPU總時(shí)間的7%，10%，16%，19%，所以可以得出結(jié)論：計(jì)算基因型概率時(shí)間約為總時(shí)間的20%。

通過(guò)上述數(shù)據(jù)可以發(fā)現(xiàn)ANGSD軟件擁有進(jìn)一步提速的可能性，接下來(lái)本文將對(duì)其軟件加速的方案進(jìn)行探討分析。

3 兩種并行方式簡(jiǎn)介

3.1 OpenMP并行

并行計(jì)算機(jī)可簡(jiǎn)單分成共享內(nèi)存和分布式存儲(chǔ)器。共享內(nèi)存是指幾個(gè)內(nèi)核之間共用的一條存儲(chǔ)器。而對(duì)只注重計(jì)算、只需要對(duì)線(xiàn)程間關(guān)系做出最基本控制的應(yīng)用而言，OpenMP十分適用。OpenMP屬于半自動(dòng)多核并行，它介于手工并行和全自動(dòng)并行之間。它通過(guò)兩方面實(shí)現(xiàn)并行化：（1）依靠編譯器來(lái)實(shí)現(xiàn)多核并行化；（2）需要程序員指定并行的代碼段或分配的粒度。

通常情況下，人工介入的部分只需要加上注釋或是標(biāo)注，相比手工并行，半自動(dòng)并行更加通俗易懂［3］ 。

OpenMP是一個(gè)應(yīng)用程序接口（API），由一組主要的計(jì)算機(jī)硬件和軟件供應(yīng)商聯(lián)合定義。OpenMP為共享內(nèi)存并行應(yīng)用程序的開(kāi)發(fā)人員提供了一個(gè)可移植、可伸縮的模型。該API在多種體系結(jié)構(gòu)上支持C/C++和Fortran。

OpenMP還給出了對(duì)并行算法的高層抽象說(shuō)明，更適用于在多核CPU計(jì)算機(jī)上的并行程序?qū)崿F(xiàn)。它在多核心CPU的機(jī)器中采用OpenMP程序方式實(shí)現(xiàn)了并行計(jì)算，編譯時(shí)通過(guò)在程序中附加的pragma命令，可以自動(dòng)地對(duì)程序并行處理，使OpenMP減少了對(duì)并行程序?qū)崿F(xiàn)的困難度和復(fù)雜性［4］。

3.2 Pthreads并行

POSIX線(xiàn)程（POSIX Threads，常被縮寫(xiě)為Pthr-eads）是POSIX的線(xiàn)程標(biāo)準(zhǔn)，定義了創(chuàng)建和操縱線(xiàn)程的一套API。Pthreads作為手工多核并行，程序員需要在程序中手動(dòng)調(diào)用API來(lái)實(shí)現(xiàn)并行化。

學(xué)術(shù)界將符合POSIX 線(xiàn)程標(biāo)準(zhǔn)的庫(kù)稱(chēng)為Pthreads，一般應(yīng)用于Unix-likePOSIX 系統(tǒng)，例如Linux系統(tǒng)。但是Windows上的實(shí)現(xiàn)也存在，例如直接使用Windows API實(shí)現(xiàn)的第三方庫(kù)pthreads-w32；而利用Windows的SFU/SUA子系統(tǒng)，則可以使用微軟提供的一部分原生POSIX API 。

Pthreads API中大致共有100個(gè)函數(shù)調(diào)用，全都以“pthread_”開(kāi)頭，并可以分為4類(lèi)：

（1）線(xiàn)程管理，例如創(chuàng)建線(xiàn)程、等待（join）線(xiàn)程、查詢(xún)線(xiàn)程狀態(tài)等。

（2）互斥鎖（Mutex）：創(chuàng)建、摧毀、鎖定、解鎖、設(shè)置屬性等操作。

（3）條件變量（Condition Variable）：創(chuàng)建、摧毀、等待、通知、設(shè)置與查詢(xún)屬性等操作。

（4）使用了互斥鎖的線(xiàn)程間的同步管理。

不但如此，Pthreads還擁有一套貼近底層的API，可以用來(lái)實(shí)現(xiàn)對(duì)線(xiàn)程的精準(zhǔn)控制，例如創(chuàng)建/取消/同步某一個(gè)線(xiàn)程，設(shè)置各個(gè)線(xiàn)程的屬性，操作互斥鎖和條件變量，進(jìn)行線(xiàn)程間同步等功能。

3.3 OpenMP和Pthreads共用模型

Pthreads和OpenMP都是多線(xiàn)程編程模型，這兩種并行語(yǔ)言通常采用如圖2所示的fork-join執(zhí)行。fork遞歸地將任務(wù)分解為較小的獨(dú)立子任務(wù)，直到它們足夠簡(jiǎn)單以便異步執(zhí)行。join將所有子任務(wù)的結(jié)果遞歸地連接成單個(gè)結(jié)果，或者在返回void的任務(wù)的情況下，程序只是等待每個(gè)子任務(wù)執(zhí)行完畢。

4 ANGSD軟件兩種并行方式區(qū)別

ANGSD作為高通量測(cè)序的一個(gè)軟件在基因測(cè)序領(lǐng)域有著不可忽視的作用。該程序使用的是Pthreads并行，而本文進(jìn)行對(duì)比的是OpenMP并行，雖然都是針對(duì)共享內(nèi)存而編寫(xiě)的API，但是這兩者存在著本質(zhì)上的區(qū)別。Pthreads要求程序員直接表明對(duì)一個(gè)線(xiàn)程的執(zhí)行情況；而OpenMP則相對(duì)而言要求寬松很多，它允許程序員只需簡(jiǎn)短地說(shuō)明哪條程序需要并行執(zhí)行，屆時(shí)將由編譯器和所執(zhí)行時(shí)的操作系統(tǒng)來(lái)決定具體的線(xiàn)程去執(zhí)行某個(gè)任務(wù)。這種差異也證明了為什么共享內(nèi)存程序同時(shí)具有兩種完全不同的API。不過(guò)在Pthreads所有線(xiàn)程的行為中的細(xì)節(jié)都必須由程序員確定。與此相反，OpenMP最主要的優(yōu)勢(shì)之一，便是它的設(shè)計(jì)使編寫(xiě)者能夠逐步地并行現(xiàn)有的串聯(lián)型程式，而不會(huì)再?gòu)牧汩_(kāi)始寫(xiě)并行程序。這可以極大地簡(jiǎn)化編程煩瑣的程序，因此把Pthreads更換為OpenMP并行進(jìn)行簡(jiǎn)化是可行的［5］ 。

5 結(jié)語(yǔ)

本文探討的內(nèi)容對(duì)ANGSD軟件加速的可能性進(jìn)行了初步的分析，ANGSD軟件執(zhí)行的流程圖以及Bam文件在虛擬機(jī)以及服務(wù)器的運(yùn)行也證明了其軟件是有加速的可行性的。對(duì)比Pthreads和OpenMP的優(yōu)缺點(diǎn)也可以說(shuō)明OpenMP相對(duì)Pthreads而言具有更為簡(jiǎn)化的編程方式，因此接下來(lái)的工作就可以對(duì)ANGSD軟件進(jìn)行程序的更改，提高其運(yùn)行的速度。

參考文獻(xiàn)

［1］黃芝準(zhǔn).組學(xué)大數(shù)據(jù)環(huán)境下的基因信息并行處理與分析方法研究［D］.合肥：中國(guó)科學(xué)技術(shù)大學(xué)，2017.

［2］DURVASULA A，HOFFMAN P J，KENT T V，et al. Angsd-wrapper：utilities for analysing next-generation sequencing data.［J］. Molecular ecology resources，2016（6）：1449-1450.

［3］戴晨，陳鵬，楊冬蕾，等.面向多核的并行編程和優(yōu)化研究［J］.計(jì)算機(jī)應(yīng)用與軟件，2013（12）：198-202，279.

［4］產(chǎn)院東.基于多核和眾核平臺(tái)的并行DNA序列比對(duì)算法［D］.濟(jì)南：山東大學(xué)，2019.

［5］段皞一.MPI、OpenMP、Taichi并行編程語(yǔ)言探究［J］.電子元器件與信息技術(shù)，2022（4）：123-134.

（編輯 王永超）

Abstract： In recent years， multi-core technology has been widely used in gene sequencing， so the parallelization of gene sequencing has become a hot topic in China. In the field of gene sequencing， with the rapid development and application of high-throughput sequencing technology， a large amount of data is generated all the time， and the single program operation of this serial program cannot meet the needs of society for gene sequencing. In order to analyze thousands of samples simultaneously， we need fast， flexible and memory-efficient implementations. OpenMP program parallelism for a program called ANGSD is very necessary.

Key words： OpenMP; next-generation sequencing; multithreading