本發(fā)明涉及生物和計算機領域，特別涉及通過大數(shù)據(jù)技術高通量測序數(shù)據(jù)統(tǒng)計方法和裝置。

背景技術：

基于高通量測序數(shù)據(jù)尋找和致病基因、癌癥治療、個性化用藥相關的染色體突變位點為臨床應用提供了不可估量的前景。由于測序技術的不斷進步，獲取到的數(shù)據(jù)越來越多，如何快速地處理不斷增加的高通量測序數(shù)據(jù)已成為亟待解決的問題。

在獲取到高通量測序數(shù)據(jù)后，需要對原始數(shù)據(jù)進行統(tǒng)計。目前常用的軟件工具是FastQC，但是FastQC只能在單臺機器上運行，運行速度較慢、運行時間較長，如一個3.8G的50基因的高通量測序數(shù)據(jù)在FastQC上需要運行6分鐘以上。隨著數(shù)據(jù)量的增加，F(xiàn)astQC消耗的處理的時間也越來越長。因此急需縮短高通量測序數(shù)據(jù)在統(tǒng)計環(huán)節(jié)消耗的時間，使高質量的數(shù)據(jù)能夠快速地進入后續(xù)分析流程。

技術實現(xiàn)要素：

有鑒于此，本發(fā)明基于分布式計算框架提供了一種對高通量測序數(shù)據(jù)速度更快的統(tǒng)計方法和統(tǒng)計裝置。

本發(fā)明的實施例提供了一種對高通量測序數(shù)據(jù)的統(tǒng)計方法，所述方法包括：

根據(jù)高通量測序數(shù)據(jù)為并行計算做準備；

對準備好的高通量測序數(shù)據(jù)進行并行計算；

將并行計算結果匯總得出統(tǒng)計數(shù)據(jù)。

優(yōu)選地，所述根據(jù)所述高通量測序數(shù)據(jù)為并行計算做準備包括：

根據(jù)所述高通量測序數(shù)據(jù)中的堿基質量值確定堿基質量值轉換方式；

對已確定堿基質量值轉換方式的高通量測序數(shù)據(jù)進行切分；

生成對切分后的數(shù)據(jù)塊進行并行計算的執(zhí)行實體。

優(yōu)選地，所述對已確定堿基質量值轉換方式的高通量測序數(shù)據(jù)進行切分包括：

將包含已確定堿基質量值轉換方式的文件轉換為RDD；

將RDD切分為partition。

優(yōu)選地，所述生成對切分后的數(shù)據(jù)塊進行并行計算的執(zhí)行實體包括：生成對partition進行并行計算的執(zhí)行實體task。

優(yōu)選地，所述對準備好的高通量測序數(shù)據(jù)進行并行計算包括：通過執(zhí)行實體并行地計算每個切分后的數(shù)據(jù)塊中與序列行相關的統(tǒng)計信息以及堿基質量值。

優(yōu)選地，所述將并行計算結果匯總得出統(tǒng)計數(shù)據(jù)包括：根據(jù)每個切分后的數(shù)據(jù)塊中與序列行相關的統(tǒng)計信息和堿基質量值統(tǒng)計所述高通量測序數(shù)據(jù)中每一列的堿基質量值分布。

本發(fā)明的實施例還提供了一種對高通量測序數(shù)據(jù)的統(tǒng)計裝置，所述裝置包括：

并行準備模塊，用于根據(jù)高通量測序數(shù)據(jù)為并行計算做準備；

并行計算模塊，用于對準備好的高通量測序數(shù)據(jù)進行并行計算；

結果匯總模塊，用于將并行計算結果匯總得出統(tǒng)計數(shù)據(jù)。

優(yōu)選地，所述并行準備模塊包括：

堿基質量值轉換方式確定單元，用于根據(jù)所述高通量測序數(shù)據(jù)中的堿基質量值確定堿基質量值轉換方式；

數(shù)據(jù)切分單元：用于對已確定堿基質量值轉換方式的高通量測序數(shù)據(jù)進行切分；

執(zhí)行實體生成單元：用于生成對切分后的數(shù)據(jù)塊進行并行計算的執(zhí)行實體。

優(yōu)選地，所述數(shù)據(jù)切分單元具體用于：

將包含已確定堿基質量值轉換方式的文件轉換為RDD；

將RDD切分為partition。

優(yōu)選地，所述執(zhí)行實體生成單元具體用于：生成對partition進行并行計算的執(zhí)行實體task。

優(yōu)選地，所述并行計算模塊具體用于：通過執(zhí)行實體并行地計算每個切分后的數(shù)據(jù)塊中與序列行相關的統(tǒng)計信息以及堿基質量值。

優(yōu)選地，所述結果匯總模塊具體用于：根據(jù)每個切分后的數(shù)據(jù)塊中與序列行相關的統(tǒng)計信息和堿基質量值統(tǒng)計所述高通量測序數(shù)據(jù)中每一列的堿基質量值分布。

本發(fā)明提供的技術方案采用大數(shù)據(jù)處理技術Spark集群開發(fā)了用于對高通量測序數(shù)據(jù)進行統(tǒng)計的軟件工具SfastQC，極大地提高了對高通量測序數(shù)據(jù)進行統(tǒng)計的速度：在采用12核CPU、128G內存的單機虛擬三個節(jié)點的條件下，SfastQC處理3.8G的50基因的測序數(shù)據(jù)所需時間縮短為42秒，比目前的FastQC快8倍以上。如果機器配置更高，數(shù)據(jù)量更大，那么加速效果會更加明顯。

附圖說明

圖1為本發(fā)明實施例一提供的高通量測序數(shù)據(jù)統(tǒng)計方法的示意圖；

圖2為本發(fā)明實施例二提供的高通量測序數(shù)據(jù)統(tǒng)計方法的示意圖；

圖3為本發(fā)明實施例三提供的高通量測序數(shù)據(jù)統(tǒng)計方法的示意圖；

圖4為本發(fā)明實施例三提供的高通量測序數(shù)據(jù)統(tǒng)計方法繪制的箱線圖；

圖5為本發(fā)明實施例四提供的高通量測序數(shù)據(jù)統(tǒng)計裝置的示意圖；

圖6為本發(fā)明實施例五、六提供的高通量測序數(shù)據(jù)統(tǒng)計方法的示意圖；

具體實施方式

為使本發(fā)明的目的、技術方案和優(yōu)點更加清楚，以下將參照本發(fā)明實施例中的附圖，通過實施方式清楚、完整地描述本發(fā)明的技術方案，顯然，所描述的實施例是本發(fā)明一部分實施例，而不是全部的實施例。

在本發(fā)明的實施例中，相關名詞解釋如下：

高通量測序數(shù)據(jù)：用高通量測序方法獲得的數(shù)據(jù)；

高通量測序：能夠一次并行對幾十萬到幾百萬條DNA片段進行序列測定的測序方法；

第一閾值：在本發(fā)明的實施例中為58，在實際應用中可以根據(jù)具體情況取適當?shù)闹担?/p>

第二閾值：在本發(fā)明的實施例中為76，在實際應用中可以根據(jù)具體情況取適當?shù)闹担?/p>

第一堿基質量值轉換方式：將質量行中每一位置上的ASCII碼字符的ASCII碼值減去第一閾值，得到對應堿基的質量值；

第二堿基質量值轉換方式：將質量行中每一位置上的ASCII碼字符的ASCII碼值減去第二閾值，得到對應堿基的質量值。

實施例一

請參閱圖1，在本發(fā)明的第一個實施例中，在根據(jù)獲得的高通量測序數(shù)據(jù)為并行計算做好準備后，對數(shù)據(jù)進行并行計算，根據(jù)計算結果繪制統(tǒng)計圖。

S101、根據(jù)高通量測序數(shù)據(jù)為并行計算做準備。

在包含高通量測序數(shù)據(jù)的FastQ文件中，每一條記錄包括四行，分別為：

以“@”開頭后面附加測序介紹信息的標識行；

由A、T、G、C四種堿基組成的序列行(測序儀無法識別的堿基用N表示)；

“+”行(或者“+”后面附帶標識行中@后面的內容，但該內容一般被省略)；

由ASCII碼字符組成的質量行(質量行和序列行長度相等且質量行中的ASCII碼字符與序列行中的堿基一一對應，質量行中每一位置的ASCII碼字符代表與其對應的序列行中對應位置的堿基質量值)。

將包含原始高通量測序數(shù)據(jù)的FastQ文件切分為多個數(shù)據(jù)塊，為并行計算做準備。

S102、對準備好的高通量測序數(shù)據(jù)進行并行計算。

根據(jù)實際需求對多個數(shù)據(jù)塊進行并行計算，例如：計算每個堿基的質量值、高通量測序數(shù)據(jù)的序列行長度等等。

S103、將并行計算結果匯總得出統(tǒng)計數(shù)據(jù)。

將多個數(shù)據(jù)塊并行計算得出的多個中間結果匯總，得到最終統(tǒng)計數(shù)據(jù)，并根據(jù)這些統(tǒng)計數(shù)據(jù)繪制相應的高通量測序數(shù)據(jù)統(tǒng)計圖。

實施例二

請參閱圖2，在本發(fā)明的第二個實施例中，使用Hadoop并行計算框架對高通量測序數(shù)據(jù)進行并行計算，將多臺計算機上并行計算的結果進行合并，得到最終統(tǒng)計數(shù)據(jù)并繪制統(tǒng)計圖。

在本實施例中，相關名詞解釋如下：

Hadoop：由Apache基金會所開發(fā)的分布式并行計算框架。

HDFS(Hadoop Distributed File System)：由Hadoop實現(xiàn)的一個分布式文件系統(tǒng)。

S201、根據(jù)所述高通量測序數(shù)據(jù)中的堿基質量值確定堿基質量值轉換方式。

輸入包含原始高通量測序數(shù)據(jù)的FastQ文件，在FastQ文件中，每一條記錄包括四行，其中序列行由A、T、G、C四種堿基組成(測序儀無法識別的堿基用N表示)；質量行由ASCII碼字符組成，質量行和序列行長度相等，質量行中每一位置上的ASCII碼字符代表與其對應的序列行中的堿基的質量。

根據(jù)讀入的FastQ文件質量行中的堿基質量值確定本文件對應的堿基質量值轉換方式：

如果讀取到大于第一閾值且小于等于第二閾值的堿基質量值，則忽略此值，繼續(xù)讀入下一位置的堿基質量值；

如果讀取到小于等于第一閾值的堿基質量值，則確定本文件對應第一堿基質量值轉換方式，選擇堿基質量值轉換方式的過程結束；

如果讀取到大于第二閾值的堿基質量值，則確定本文件對應第二堿基質量值轉換方式，選擇堿基質量值轉換方式的過程結束。

S202、對已確定堿基質量值轉換方式的高通量測序數(shù)據(jù)進行切分。

HDFS將輸入的FastQ文件根據(jù)一定的規(guī)則切分成小數(shù)據(jù)塊并保存。切分規(guī)則如下：例如輸入文件為3G，在Hadoop中將一個數(shù)據(jù)塊的大小設置為128M，那么輸入文件總共將被切分為3*1024/128＝24塊。

S203、生成包含切分后的數(shù)據(jù)的執(zhí)行實體map task。

在Hadoop集群中，參與并行計算的多臺計算機并行地讀取HDFS中的小數(shù)據(jù)塊，并啟動一個job，job為每一個小數(shù)據(jù)塊生成一個map task。Map task是并行計算的執(zhí)行實體。

S204、對準備好的高通量測序數(shù)據(jù)進行并行計算。

Hadoop在參與并行計算的多臺計算機上并行地運行map task。Map task以小數(shù)據(jù)塊為單位進行統(tǒng)計(如統(tǒng)計每個小數(shù)據(jù)塊中的堿基總數(shù)，或堿基的GC含量等)，并根據(jù)實際需要對每個小數(shù)據(jù)塊進行計算(如計算堿基質量值等)，然后把對每個小數(shù)據(jù)塊的統(tǒng)計和計算結果寫到HDFS中。

S205、將并行計算結果匯總得出統(tǒng)計數(shù)據(jù)。

Hadoop根據(jù)用戶指定的數(shù)量生成多個reduce task。Reduce task從HDFS中讀取map task對每個小數(shù)據(jù)塊的統(tǒng)計和計算結果，并將其進行合并，得出最終統(tǒng)計數(shù)據(jù)，然后根據(jù)最終統(tǒng)計數(shù)據(jù)繪制相應的高通量測序數(shù)據(jù)統(tǒng)計圖。

實施例三

請參閱圖3，在本發(fā)明的第三個實施例中，使用spark并行計算框架對高通量測序數(shù)據(jù)進行并行計算，將多臺計算機上并行計算的結果進行合并，得到最終統(tǒng)計數(shù)據(jù)并繪制統(tǒng)計圖。

S301、根據(jù)所述高通量測序數(shù)據(jù)中的堿基質量值確定堿基質量值轉換方式。

輸入包含原始高通量測序數(shù)據(jù)的FastQ文件，在FastQ文件中，每一條記錄包括四行，其中序列行由A、T、G、C四種堿基組成(測序儀無法識別的堿基用N表示)；質量行由ASCII碼字符組成，質量行和序列行長度相等，質量行中每一位置上的ASCII碼字符代表與其對應的序列行中對應位置的堿基質量值。

根據(jù)讀入的FastQ文件質量行中的堿基質量值確定本文件對應的堿基質量值轉換方式：

如果讀取到大于第一閾值且小于等于第二閾值的堿基質量值，則忽略此值，繼續(xù)讀入下一位置的堿基質量值；

如果讀取到小于等于第一閾值的堿基質量值，則確定本文件對應第一堿基質量值轉換方式，選擇堿基質量值轉換方式的過程結束；

如果讀取到大于第二閾值的堿基質量值，則確定本文件對應第二堿基質量值轉換方式，選擇堿基質量值轉換方式的過程結束。

S302、將包含已確定堿基質量值轉換方式的文件轉換為RDD；將RDD切分為partition。

在讀取FastQ文件時，spark先將其轉換為RDD，在生成RDD時用戶可以根據(jù)實際需要指定將RDD切分為partition的數(shù)量。例如輸入文件為3G，設置將RDD切分為24個partition，那么每一個partition所占存儲空間為3*1024/24＝128M。實際生成的partition的數(shù)量最少為(該文件所占存儲空間/128M)，如果指定的partition數(shù)量少于(該文件所占存儲空間/128M)，則實際將生成(該文件所占存儲空間/128M)個partition。Spark根據(jù)用戶指定的partition數(shù)量將RDD切分成若干partition。

在本實施例中：

Spark：是UC Berkeley AMPLab開發(fā)的一種計算框架。

RDD是指彈性分布式數(shù)據(jù)集(Resilient Distributed Datasets)，它是可容錯的并行數(shù)據(jù)結構，使用戶能夠顯式地在內存中保存中間的運算結果，通過控制RDD的分區(qū)來優(yōu)化數(shù)據(jù)的布局，并使用豐富的轉換算子進行操作。

Partition是指spark在計算過程中，生成的數(shù)據(jù)在計算空間內的最小單元。

S303、生成對partition進行并行計算的執(zhí)行實體task。

在本實施例中：

Job是指包在spark中含由多個stage組成的并行計算，對RDD執(zhí)行action操作后會生成job；

Stage是指在spark中，一個job會根據(jù)處理過程的需要而分成不同的階段即stage，stage由多個task組成；

Task是指被送到為某個應用啟動的executor進程的工作單元。

在spark中對RDD進行action操作時生成DAG Scheduler(有向無環(huán)圖調度器)，從而啟動一個job。對一個job內的操作，根據(jù)處理過程是否需要shuffle分成不同的stage，并在每一個stage內產(chǎn)生一系列的task。通常一個RDD內的task數(shù)量與partition的數(shù)量相同。后續(xù)多個執(zhí)行實體task將在多臺計算機上對不同的partition執(zhí)行并行計算過程。

S304、通過執(zhí)行實體task并行地計算每個partition中與序列行相關的統(tǒng)計信息以及堿基質量值。

在多臺計算機上通過多個執(zhí)行實體task并行地對每個partition進行以下計算：

對高通量測序數(shù)據(jù)中每個記錄的序列行計數(shù)得到該記錄的序列總數(shù)、計算每個記錄的序列行中的堿基總數(shù)和堿基中的GC含量；

計算堿基質量值：

如果本次并行計算數(shù)據(jù)對應的是第一堿基質量值轉換方式，則將質量行中每一位置上的ASCII碼字符的ASCII碼值減去第一閾值即為對應的堿基質量值；

如果本次并行計算數(shù)據(jù)對應的是第二堿基質量值轉換方式，則將質量行中每一位置上的ASCII碼字符的ASCII碼值減去第二閾值，即為對應的堿基質量值。

S305、根據(jù)每個partition中與序列行相關的統(tǒng)計信息和堿基質量值統(tǒng)計所述高通量測序數(shù)據(jù)中每一列的堿基質量值分布。

將所有partition中的序列數(shù)累加得出序列總數(shù)；統(tǒng)計所有partition中的序列長度得出序列長度范圍；將所有partition中的堿基數(shù)累加得出堿基總數(shù)；將所有partition中的堿基質量值累加并除以堿基總數(shù)得出每一列堿基質量平均值；將所有partition中的堿基G和C的數(shù)量累加并除以堿基總數(shù)得到堿基GC含量。

根據(jù)每個partition中的計算結果統(tǒng)計高通量測序數(shù)據(jù)中每一列的堿基質量值分布：10％點、上四分位數(shù)(25％點)、中位數(shù)(50％點)、下四分位數(shù)(75％點)和90％點。

根據(jù)上述統(tǒng)計數(shù)據(jù)生成結果文件，此結果文件包括原始高通量測序數(shù)據(jù)文件的名稱、序列總數(shù)、序列長度范圍、堿基總數(shù)、平均堿基質量值、GC堿基含量和每一列的堿基質量值分布；根據(jù)每一列的堿基質量值分布即每一列堿基質量值的10％點、上四分位數(shù)(25％點)、中位數(shù)(50％)、下四分位數(shù)(75％點)、90％點和平均值繪制箱線圖(箱線圖是利用數(shù)據(jù)中的上述五個統(tǒng)計量來描述數(shù)據(jù)的一種方法，根據(jù)箱線圖可以大致看出數(shù)據(jù)是否具有有對稱性，分布的分散程度等信息，特別可以用于對幾個樣本的比較)。如圖4所示例，在本實施例的箱線圖中，縱坐標為堿基質量值，分布在0到42之間；橫坐標為堿基位置坐標，即高通量測序數(shù)據(jù)中的序列長度。

實施例四

如圖5所示，本發(fā)明的第四個實施例提供了一種高通量測序數(shù)據(jù)的統(tǒng)計裝置，所述裝置包括：

并行準備模塊510，用于根據(jù)高通量測序數(shù)據(jù)為并行計算做準備；

并行計算模塊520，用于對準備好的高通量測序數(shù)據(jù)進行并行計算；

結果匯總模塊530，用于將并行計算結果匯總得出統(tǒng)計數(shù)據(jù)。

在本實施例中，并行準備模塊510將包含原始高通量測序數(shù)據(jù)的FastQ文件切分為多個數(shù)據(jù)塊，為并行計算做準備。并行計算模塊520根據(jù)實際需求對多個數(shù)據(jù)塊進行并行計算。將多個數(shù)據(jù)塊并行計算得出的多個中間結果匯總，得到最終統(tǒng)計數(shù)據(jù)，結果匯總模塊530并根據(jù)這些統(tǒng)計數(shù)據(jù)繪制相應的高通量測序數(shù)據(jù)統(tǒng)計圖。

實施例五

如圖6所示，本發(fā)明的第五個實施例提供了一種高通量測序數(shù)據(jù)的統(tǒng)計裝置，所述裝置包括：

并行準備模塊610，用于根據(jù)高通量測序數(shù)據(jù)為并行計算做準備；

并行計算模塊620，具體用于通過執(zhí)行實體并行地計算每個切分后的數(shù)據(jù)塊中與序列行相關的統(tǒng)計信息以及堿基質量值。

結果匯總模塊630，具體用于根據(jù)每個切分后的數(shù)據(jù)塊中與序列行相關的統(tǒng)計信息和堿基質量值統(tǒng)計所述高通量測序數(shù)據(jù)中每一列的堿基質量值分布。

所述并行準備模塊610包括：

堿基質量值轉換方式確定單元6101，用于根據(jù)所述高通量測序數(shù)據(jù)中的堿基質量值確定堿基質量值轉換方式；

數(shù)據(jù)切分單元6102：用于對已確定堿基質量值轉換方式的高通量測序數(shù)據(jù)進行切分；

執(zhí)行實體生成單元6103：用于生成對切分后的數(shù)據(jù)塊進行并行計算的執(zhí)行實體。

在本實施例中，并行準備模塊610中的堿基質量值轉換方式確定單元6101根據(jù)輸入的高通量測序數(shù)據(jù)中的堿基質量值確定堿基質量值轉換方式。并行準備模塊610的數(shù)據(jù)切分單元6102對已確定堿基質量值轉換方式的高通量測序數(shù)據(jù)進行切分。并行準備模塊610中的執(zhí)行實體生成單元6103生成包含切分后的數(shù)據(jù)的執(zhí)行實體map task。并行計算模塊620通過map task以小數(shù)據(jù)塊為單位并行地進行統(tǒng)計和計算。結果匯總模塊630通過reduce task將多個計算機上對每個小數(shù)據(jù)塊的統(tǒng)計和計算結果進行合并，得到最終統(tǒng)計數(shù)據(jù)，并根據(jù)最終統(tǒng)計數(shù)據(jù)繪制相應的高通量測序數(shù)據(jù)統(tǒng)計圖。

實施例六

如圖6所示，本發(fā)明的第六個實施例提供了一種高通量測序數(shù)據(jù)的統(tǒng)計裝置，所述裝置包括：

并行準備模塊610，用于根據(jù)高通量測序數(shù)據(jù)為并行計算做準備；

并行計算模塊620，用于對準備好的高通量測序數(shù)據(jù)進行并行計算；

結果匯總模塊630，用于將并行計算結果匯總得出統(tǒng)計數(shù)據(jù)。

所述并行準備模塊610包括：

堿基質量值轉換方式確定單元6101，用于根據(jù)所述高通量測序數(shù)據(jù)中的堿基質量值確定堿基質量值轉換方式；

數(shù)據(jù)切分單元6102，具體用于將包含已確定堿基質量值轉換方式的文件轉換為RDD；將RDD切分為partition。

執(zhí)行實體生成單元6103，用于生成對partition進行并行計算的執(zhí)行實體。

所述并行計算模塊620，具體用于通過執(zhí)行實體task并行地計算每個partition中與序列行相關的統(tǒng)計信息以及堿基質量值。

所述結果匯總模塊630，具體用于根據(jù)每個partition中與序列行相關的統(tǒng)計信息和堿基質量值統(tǒng)計所述高通量測序數(shù)據(jù)中每一列的堿基質量值分布。

在本發(fā)明的第六個實施例中，并行準備模塊610中的堿基質量值轉換方式確定單元6101根據(jù)輸入的高通量測序數(shù)據(jù)中的堿基質量值確定堿基質量值轉換方式。并行準備模塊610的數(shù)據(jù)切分單元6102將包含已確定堿基質量值轉換方式的文件轉換為RDD；將RDD切分為partition。并行準備模塊610中的執(zhí)行實體生成單元6103生成對partition進行并行計算的執(zhí)行實體task。并行計算模塊620通過執(zhí)行實體task并行地計算每個partition中與序列行相關的統(tǒng)計信息以及堿基質量值。結果匯總模塊630根據(jù)每個partition中與序列行相關的統(tǒng)計信息和堿基質量值統(tǒng)計所述高通量測序數(shù)據(jù)中每一列的堿基質量值分布，并繪制箱線圖。

本領域普通技術人員可以理解，實現(xiàn)上述本發(fā)明實施例中的高通量測序數(shù)據(jù)統(tǒng)計方法和統(tǒng)計裝置可以通過程序指令相關的硬件來完成，所述的程序可以存儲于可讀取存儲介質中，該程序在執(zhí)行時執(zhí)行上述方法中的對應步驟。

以上所述僅是本發(fā)明的優(yōu)選實施方式，應當指出，對于本技術領域的普通技術人員來說，在不脫離本發(fā)明原來的前提下，還可以做出若干改進和潤飾，這些改進和潤飾也應視為本發(fā)明的保護范圍。