本文關(guān)鍵詞：基因組二代測序數(shù)據(jù)的自動化分析流程，由筆耕文化傳播整理發(fā)布。

Hereditas (Beijing) 2014年6月, 36(6): 618―624

技術(shù)與方法

基因組二代測序數(shù)據(jù)的自動化分析流程

李文軻1, 李豐余1,2, 張思瑤1, 蔡斌1, 鄭娜1, 聶宇1, 周到2, 趙倩1

1. 中國醫(yī)學(xué)科學(xué)院, 北京協(xié)和醫(yī)學(xué)院, 國家心血管病中心, 阜外心血管病醫(yī)院, 心血管疾病國家重點實驗室, 北京 100037; 2. 中南民族大學(xué)生物醫(yī)學(xué)工程學(xué)院, 武漢430074

摘要： 二代測序技術(shù)的發(fā)展對測序數(shù)據(jù)的處理分析提出了很高的要求。目前二代測序數(shù)據(jù)分析軟件很多, 但是

絕大多數(shù)軟件僅能完成單一的分析功能(例如：僅進行序列比對或變異讀取或功能注釋等), 如何能正確高效地選擇整合這些軟件已成為迫切需求。文章設(shè)計了一套基于perl語言和SGE資源管理的自動化處理流程來分析Illumina平臺基因組測序數(shù)據(jù)。該流程以測序原始序列數(shù)據(jù)作為輸入, 調(diào)用業(yè)界標(biāo)準(zhǔn)的數(shù)據(jù)處理軟件(如：BWA, Samtools, GATK, ANNOVAR等), 最終生成帶有相應(yīng)功能注釋、便于研究者進一步分析的變異位點列表。該流程通過自動化并行腳本控制流程的高效運行, 一站式輸出分析結(jié)果和報告, 簡化了數(shù)據(jù)分析過程中的人工操作, 大大提高了運行效率。用戶只需填寫配置文件或使用圖形界面輸入即可完成全部操作。該工作為廣大研究者分析二代測序數(shù)據(jù)提供了便利的途徑。

關(guān)鍵詞：　二代測序; 自動化數(shù)據(jù)分析; 流程; 變異檢測

Automatic analysis pipeline of next-generation sequencing data

Wenke Li1, Fengyu Li1, 2, Siyao Zhang1, Bin Cai1, Na Zheng1, Yu Nie1, Dao Zhou2, Qian Zhao1

1. State Key Laboratory of Cardiovascular Disease, Fuwai Hospital, National Center for Cardiovascular Disease, Chinese Academy of Medical Sciences and Peking Union Medical College, Beijing 100037, China;

2. College of Biomedical Engineering, South-Central University for Nationalities, Wuhan 430074, China

Abstract: The development of next-generation sequencing has generated high demand for data processing and analysis. Although there are a lot of software for analyzing next-generation sequencing data, most of them are designed for one specific function (e.g., alignment, variant calling or annotation). Therefore, it is necessary to combine them together for data analysis and to generate interpretable results for biologists. This study designed a pipeline to process Illumina sequencing data based on Perl programming language and SGE system. The pipeline takes original sequence data (fastq format) as input, calls the standard data processing software (e.g., BWA, Samtools, GATK, and Annovar), and finally outputs a list of annotated vari-ants that researchers can further analyze. The pipeline simplifies the manual operation and improves the efficiency by

收稿日期: 2013?09?07; 修回日期: 2014?01?20

基金項目:國家重點基礎(chǔ)研究發(fā)展計劃(973計劃)項目(編號：2010CB529505)和中央高校基本科研業(yè)務(wù)費專項資金(編號：2012-XHGX02)資助

作者簡介:李文軻, 碩士, 助理研究員, 研究方向：生物信息學(xué)。Tel: 010-88396071; E-mail: wksofia@gmail.com

通訊作者:趙倩, 博士, 副研究員, 研究方向：遺傳學(xué), 生物信息學(xué)。E-mail: zhaoqian82@gmail.com DOI: 10.3724/SP.J.1005.2014.0618 網(wǎng)絡(luò)出版時間: 2014-3-25 17:22:30

URL:

第6期

李文軻等: 基因組二代測序數(shù)據(jù)分析流程的自動化實現(xiàn)

619

automatization and parallel computation. Users can easily run the pipeline by editing the configuration file or clicking the graphical interface. Our work will facilitate the research projects using the sequencing technology.

Keywords: next generation sequencing; automatic data analysis; pipeline; variantion detection

二代測序技術(shù)(Next-generation sequencing)大幅度降低了測序的時間和成本, 使得大規(guī)模測序逐漸成為常規(guī)的實驗室研究和臨床檢測手段。測序產(chǎn)生的數(shù)據(jù)量急劇增加, 如何高效地分析這些數(shù)據(jù), 已成為迫切需要解決的問題。目前, 分析序列信息的生物信息學(xué)軟件紛繁復(fù)雜, 但基本上每個軟件只能完成單一的分析功能, 實現(xiàn)一個完整的分析流程則需要對眾多軟件進行整合, 而手動串聯(lián)的效率往往不盡人意; 同時, 這些軟件需要在Linux工作環(huán)境下以命令行運行, 要求用戶具備較好的計算機背景; 另外, 即便一些實驗室完成了分析流程的構(gòu)建, 他們往往不會公開許多細節(jié), 新用戶仍然要從頭建起。本研究致力于構(gòu)建經(jīng)典的二代測序數(shù)據(jù)分析流程, 并實現(xiàn)各個環(huán)節(jié)的高效自動化管理和分析, 減輕研究者前期的工作負擔(dān), 促進相關(guān)領(lǐng)域進一步對基因組測序研究項目的順利開展。

(Hiseq2500), 甚至更高的250個堿基(Miseq)。測序讀長不斷增加, 測序通量也在不斷上升。Illumina Hiseq2500是目前世界上通量最高的測序平臺, 最多可以在大約10 d的時間內(nèi)測定3000億個堿基——即6～7個人類全基因組或60～80個人類全外顯子組的序列測定。

Illumina平臺以FASTQ格式[2]存儲測序結(jié)果, 這也是本流程的輸入文件。FASTQ文件記錄內(nèi)容包括所測的堿基讀段和質(zhì)量, 其數(shù)據(jù)格式如圖1所示。每條讀段(reads)占四行：第一行和第三行為讀段識別碼, 包含測序儀SN號、產(chǎn)生讀段的巷道(lane)、該讀段的編號等信息; 第二行為讀段測到的堿基序列; 第四行為所測到堿基的質(zhì)量分數(shù), 每一個堿基都會對應(yīng)一個質(zhì)量分數(shù)。

1.2 數(shù)據(jù)處理流程及軟件簡介

目前測序數(shù)據(jù)處理軟件很多, 我們綜合考慮了適用性和效率, 整合出了一套標(biāo)準(zhǔn)的數(shù)據(jù)處理流程。具體來說, 獲得FASTQ格式的原始測序數(shù)據(jù)后, 需要對數(shù)據(jù)進行以下處理：(1)使用BWA軟件把這些短序列和參考基因組進行對比, 確定短序列在基因組上的位置, 把短序列組裝成完整的人類參考基因組; (2)使用Samtools軟件把這些短序列調(diào)整成按一定順序(1~22, X, Y, 其他)排列的序列, 并進行數(shù)據(jù)格式的轉(zhuǎn)換; (3)使用Picard軟件把測序產(chǎn)生的冗

1 數(shù)據(jù)的獲取和分析流程的構(gòu)建

1.1 Illumina測序數(shù)據(jù)

本流程適用于Illumina測序平臺產(chǎn)出的雙端(Paired ends)測序數(shù)據(jù)。Illumina測序技術(shù)采用邊合成邊測序(Sequencing by synthesis, SBS)的方法, 早期的GA測序儀測序讀長只有100個堿基, 隨著技術(shù)的改進, 目前的讀長已經(jīng)增加到150個堿基

圖1 FASTQ格式示例

620

Hereditas (Beijing) 2014

第36卷

余信息和噪聲去掉; (4)使用GATK尋找樣本測序數(shù)據(jù)與參考基因組的差異, 列出這些差異點; (5)使用Annovar對這些變異位點進行功能注釋, 得到一個易于理解的變異位點列表。處理流程圖如圖 2 所示。

的序列比對軟件, 能高效地比對短序列和參考基因組, 找到短序列在參考基因組上的位置, 該軟件最長支持至1 Mb的短序列比對。BWT方法通過B-W轉(zhuǎn)換將基因組序列按一定規(guī)則壓縮并建立索引, 再通過查找和回溯來定位讀段, 在查找時可通過堿基替代來實現(xiàn)允許的錯配。采用Burrows-Wheeler轉(zhuǎn)換的代表軟件是Bowtie和BWA。比對結(jié)果如圖3所示：界面上方是測到的短序列, 下方是短序列所比對到的參考基因組。

1.2.2 SAM文件處理軟件Samtools

讀段定位到基因組后推薦采用SAM(Sequence Alignment/Map)格式或其二進制版本BAM格式來存儲。二進制版本可大大節(jié)省存儲空間, 但不能直接用普通文本編輯工具顯示。

SAM文件處理軟件Samtools可以很好的對SAM/BAM格式數(shù)據(jù)進行操作, 因此, 本文使用它來進行數(shù)據(jù)格式轉(zhuǎn)換和排序。

1.2.3 測序噪聲去除和測序數(shù)據(jù)評價軟件Picard

對組裝好的全基因組數(shù)據(jù), 需要將過度重復(fù)測

到的數(shù)據(jù)進行剔除, 并且需要對數(shù)據(jù)質(zhì)量進行評價, Picard軟件可以很好地完成這兩項工作。 1.2.4 變異檢測軟件GATK

GATK主要用于在測序數(shù)據(jù)中尋找變異, 包括單堿基變異(SNV)、短插入缺失(INDEL), 是當(dāng)前業(yè)界用來尋找變異的主流軟件。

圖2 處理流程圖

1.2.1 讀句比對軟件BWA

BWA(Burrows-Wheeler Alignment tool)是基于Burrows-Wheeler變換(Burrows-Wheeler Transform)

圖3 比對結(jié)果示例

利用Broad Institute的IGV(Integrated Genomics Viewer)對數(shù)據(jù)進行可視化, 圖4同。

第6期

李文軻等: 基因組二代測序數(shù)據(jù)分析流程的自動化實現(xiàn)

621

變異指測序序列和參考序列的差異。如圖4所示, 參考序列上的堿基是胸腺嘧啶(T), 而測序數(shù)據(jù)上的堿基是鳥嘌呤(G), 說明此處有一個T→G 的 突變。

1.2.5 變異注釋軟件ANNOVAR

ANNOVAR是一個用于高效注釋變異的工具。注釋信息包括變異所在的染色體, 開始位置, 結(jié)束位置, 參考序列信息和觀察到的序列信息的列表。一個變異經(jīng)過ANNOVAR注釋之后, 其功能一目了然, 便于進一步的生物學(xué)分析。

務(wù)正在進行時, Perl對它進行監(jiān)控; 當(dāng)計算完成, Perl去查看它的輸出的計算結(jié)果, 并把結(jié)果作為下一個計算任務(wù)的輸入, 往計算節(jié)點上投放新的計算任務(wù)。如此循環(huán), 直到流程運行完畢。

同時, 由于每次運行的樣本不同, 數(shù)據(jù)的輸入輸出位置也有差異。如果每處理一個新的樣本, 就要對流程源碼進行大量修改, 將不利于流程的使用。為此, 本流程定義了一個配置文件(config file)。通過配置文件可以指定：流程處理的樣品名、數(shù)據(jù)輸入輸出路徑、參考序列文件, 甚至流程中涉及到的軟件的位置、軟件的運行方式; 另外, 我們還提供了對流程中主要軟件參數(shù)的修改, 以滿足高級用戶需求。每次進行一個新樣本的分析, 不需要修改主程序代碼, 只要為其創(chuàng)建一個配置文件, 主程序會自動讀取配置文件, 生成相應(yīng)的執(zhí)行代碼。

流程文件構(gòu)成如圖5所示。

2 自動化實現(xiàn)

2.1 基于Perl語言的流程設(shè)計

本數(shù)據(jù)處理流程主要使用Perl編程語言實現(xiàn)對各個軟件的高效串接和自動化操作。一項計算任

圖4 單堿基突變示例

2.2 基于資源管理軟件(SGE)的并行設(shè)計

流程的運行環(huán)境是計算機集群, 其有別于普通PC機, 一般是由一臺管理主機來協(xié)調(diào)許多計算主機

來完成大型的計算任務(wù)。根據(jù)這樣的硬件特點來設(shè)計流程, 需要考慮以下兩個問題：(1)如何讓眾多計

圖5 分析流程結(jié)構(gòu)

  本文關(guān)鍵詞：基因組二代測序數(shù)據(jù)的自動化分析流程，由筆耕文化傳播整理發(fā)布。