新一代高通量rna測序數(shù)據(jù)的處理與分析

1、本文檔下載自HYPERLINK /文檔下載網(wǎng)，內(nèi)容可能不完整，您可以點(diǎn)擊以下網(wǎng)址繼續(xù)閱讀或下載：HYPERLINK /doc/be6b0724c324a6024ee96d8b/doc/be6b0724c324a6024ee96d8b新一代高通量 RNA 測序數(shù)據(jù)的處理與分析 _新一代高通量 RNA 測序數(shù)據(jù)的處理與分析 _生物化學(xué)與生物物理進(jìn)展,37(8):834846此處圖片未下載成功此處圖片未下載成功新一代高通量RNA測序數(shù)據(jù)的處理與分析*王曦1)汪小我1)王立坤1,2)馮智星1)張學(xué)工1)*)(1)生物信息學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，清華信息科學(xué)與技術(shù)國家實(shí)驗(yàn)室(籌)生物信息學(xué)研究部，清華大學(xué)自

2、動化系，北京100084；吉林大學(xué)計算機(jī)科學(xué)與技術(shù)學(xué)院，長春130012)摘要隨著新一代高通量DNA測序技術(shù)的快速發(fā)展，RNA測序(RNA-seq)已成為基因表達(dá)和轉(zhuǎn)錄組分析新的重要手段seq技術(shù)產(chǎn)生的海量數(shù)據(jù)為生物信息學(xué)帶來了新的機(jī)遇和挑戰(zhàn)有效地對測序數(shù)據(jù)進(jìn)行針對性的生物信息學(xué)處理和分析，成為RNA-seq技術(shù)能否在科學(xué)探索中發(fā)揮重大作用的關(guān)鍵以新一代Illumina/Solexa測序平臺所產(chǎn)生的數(shù)據(jù)為例，在扼要介紹高通量RNA-seq測序流程的基礎(chǔ)上，對RNA-seq數(shù)據(jù)處理和分析的方法和現(xiàn)有軟件做一個較為全面的綜述，并對其中有待進(jìn)一步研究的問題進(jìn)行展望關(guān)鍵詞高通量RNA測序，轉(zhuǎn)錄組，基因

3、表達(dá)，數(shù)據(jù)處理與分析，生物信息學(xué)，Q6，Q7:10.3724/SP.J.1206.2010.00151學(xué)科分類號近年來，新一代高通量測序技術(shù)得到了突飛猛進(jìn)的發(fā)展，在此基礎(chǔ)上，高通量RNA測序即RNA-seq1-5也迅速發(fā)展與基因芯片技術(shù)相比，RNA-seq無需設(shè)計探針，能在全基因組范圍內(nèi)以單堿基分辨率檢測和量化轉(zhuǎn)錄片段，并能應(yīng)用于基因組圖譜尚未完成的物種6，具有信噪比高、分辨率高、應(yīng)用范圍廣等優(yōu)勢，正成為研究基因表達(dá)和轉(zhuǎn)錄組的重要實(shí)驗(yàn)手段seq為基因組學(xué)的研究帶來了高分辨率的海量數(shù)據(jù)，如何有效處理和分析這些海量數(shù)據(jù)成為這一新技術(shù)能否帶來新的科學(xué)發(fā)現(xiàn)的關(guān)鍵，一些/doc/be6b0724c324

4、a6024ee96d8b生物信息學(xué)方法與軟件也應(yīng)運(yùn)而生本文針對當(dāng)前RNA-seq應(yīng)用的現(xiàn)實(shí)情況，嘗試以Illumina/Solexa測序平臺產(chǎn)生的mRNA-seq數(shù)據(jù)為例，對RNA測序數(shù)據(jù)的產(chǎn)生過程及數(shù)據(jù)處理和分析的基本流程、關(guān)鍵方法和現(xiàn)有軟件進(jìn)行較全面的介紹，并討論RNA-seq數(shù)據(jù)分析中存在的挑戰(zhàn)格更便宜、自動化程度更高的測序技術(shù)自2005年以來，以Roche公司的454技術(shù)、Illumina公司的Solexa技術(shù)和ABI公司的SOLiD技術(shù)為標(biāo)志的新一代測序技術(shù)相繼誕生8新一代測序技術(shù)又稱作深度測序技術(shù)，主要特點(diǎn)是測序通量高、測序時間和成本顯著下降9把這種高通量測序技術(shù)應(yīng)用到由mRNA逆

5、轉(zhuǎn)錄生成的cDNA上，從而獲得來自不同基因的mRNA片段在特定樣本中的含量，這就是mRNA測序或mRNA-seq同樣原理，各種類型的轉(zhuǎn)錄本都可以用深度測序技術(shù)進(jìn)行高通量定量檢測，統(tǒng)稱作RNA-seq或RNA測序目前，在已經(jīng)推出的幾種新一代測序平臺中，Illumina/Solexa測序平臺上的RNA-seq應(yīng)用最廣，我們以此為例來綜述RNA-seq數(shù)據(jù)處理和分析的生物信息學(xué)問題和方法.高通量測序技術(shù)簡介*國家自然科學(xué)基金資助項(xiàng)目(60702002,60721003,30873464,60905013)和東南大學(xué)生物電子學(xué)國家重點(diǎn)實(shí)驗(yàn)室開放研究基金資助項(xiàng)目.*通訊聯(lián)系人./p>

6、,E-mail:zhangxg收稿日期：2010-03-25，接受日期：2010-04-30誕生于20世紀(jì)70年代的Sanger法是最早被廣泛應(yīng)用的DNA測序技術(shù)7，也是完成人類基因組計劃的基礎(chǔ)但是，由于它測序通量低，費(fèi)時費(fèi)力，科學(xué)家們一直在尋求通量更高、速度更快、價;37(8)王曦等：新一代高通量RNA測序數(shù)據(jù)的處理與分析835/Solexa測序技術(shù)的基本原理是邊合成邊測序(sequencingb/doc/be6b0724c324a6024ee96d8bysynthesis，SBS)10-12，即測序過程是以DNA單鏈為模板，在生成互補(bǔ)鏈時，利用帶熒光標(biāo)記的dNTP發(fā)出不同顏色的熒光來確定不

7、同的堿基新加入dNTP的末端被可逆的保護(hù)基團(tuán)封閉，既保證單次反應(yīng)只能加入一個堿基，又能在該堿基讀取完畢后，將保護(hù)基團(tuán)除去，使得下一個反應(yīng)可繼續(xù)進(jìn)行為了增加熒光強(qiáng)度，使之更易被成像系統(tǒng)所采集，該技術(shù)在測序之前還需要對待測片段做橋式擴(kuò)增(bridgeamplification)13(/)初期的Illumina/Solexa測序技術(shù)只能在較短的測序讀長上(2030堿基)保證較高的正確率隨著技術(shù)的改進(jìn)，目前的讀長已經(jīng)增加到100堿基以上同時，隨著雙端測序(paired-end，PE)技術(shù)的成熟，測序長度更可達(dá)到單端測序的2倍，測序通量也隨之增加這種測序技術(shù)是Solexa公司發(fā)展起來的，2007年被Il

8、lumina公司收購，因此現(xiàn)在通常被稱為Illumina/Solexa測序技術(shù)近兩年來，Illumina/Solexa測序平臺不斷升級，相繼推出了GA(GenomeAnalyzer)、GAIIx、HiSeq2000等測序儀更多關(guān)于高通量測序平臺的介紹，可以查閱相關(guān)文獻(xiàn)9,14-16RNA鄄seq測序文庫制備和測序平臺數(shù)據(jù)輸出本小節(jié)針對Illumina/Solexa測序平臺，對RNA測序文庫制備標(biāo)準(zhǔn)和平臺底層數(shù)據(jù)產(chǎn)生做一個簡單的介紹援1RNA鄄seq測序文庫制備對于mRNA-seq實(shí)驗(yàn)，從總RNA到最終的cDNA文庫制備完成主要包括以下步驟首先，用Poly(T)寡聚核苷酸從總RNA中抽取全部帶P

9、oly(A)尾的RNA，其中的主要部分就是編碼基因所轉(zhuǎn)錄的mRNA將所得RNA隨機(jī)打斷成片段，再用隨機(jī)引物和逆轉(zhuǎn)錄酶從RNA片段合成cDNA片段然后，對cDNA片段進(jìn)行末端修復(fù)并連接測序接頭(adapter)，得到將用于測序的cDNA在以上過程，將RNA隨機(jī)片段化和采用隨機(jī)引物進(jìn)行反轉(zhuǎn)錄，都是為了使所得cDNA片段較均勻地取自各個轉(zhuǎn)錄本為提高測序效率，一般還需要用電泳切膠/doc/be6b0724c324a6024ee96d8b法獲取長度范圍在200bp(依25bp)的cDNA片段，再通過RCR擴(kuò)增，得到最終的cDNA文庫在上述文庫制備過程中，如果不是只抽取帶Poly(A)尾的RNA，而是使用

10、全部的RNA，則RNA-seq測得的就是細(xì)胞中的全部轉(zhuǎn)錄本，如果把帶Poly(A)尾的RNA過濾掉，也可以得到非編碼的RNA轉(zhuǎn)錄本，如果從總RNA中只提取長度為2123個堿基左右的RNA，則得到全部的miRNA(microRNA)轉(zhuǎn)錄本,相應(yīng)的方法也稱作miRNA-seq.樣品制備最終得到的是雙鏈cDNA文庫在后續(xù)測序中，測得的每個讀段(read)隨機(jī)地來自雙鏈cDNA的某一條鏈，從讀段序列本身無法得知它是與RNA方向相同還是倒轉(zhuǎn)互補(bǔ)，在后續(xù)的讀段定位時需要兩個方向都考慮在新基因識別等應(yīng)用中，轉(zhuǎn)錄本的方向?qū)蜃⑨層葹橹匾?，需要在文庫制備和測序中保留RNA的方向信息最近有文獻(xiàn)報道了保留方向信息

11、的RNA-seq樣品制備方法17-20援2測序平臺數(shù)據(jù)輸出將RNA-seq測序文庫加入流動槽(flowcell)中的各通道(lane)，在橋式PCR擴(kuò)增后，就可以進(jìn)行測序了測序過程中，計算機(jī)軟件同步地對熒光圖像數(shù)據(jù)進(jìn)行處理，通過分析熒光信號來確定被測堿基，并給出質(zhì)量評分按照圖像上的位置坐標(biāo)，計算機(jī)程序?qū)⑼晃恢脺y得的堿基根據(jù)測序順序連成讀段(read)由于熒光圖像文件所占有的磁盤空間很大，通常GAIIx平臺一次實(shí)驗(yàn)就能產(chǎn)生上太字節(jié)(TB)的圖像文件，所以一般情況下不予保留原始的熒光圖像數(shù)據(jù)，而是只保留程序讀出的讀段數(shù)據(jù)及對應(yīng)的質(zhì)量分值，這就是多數(shù)實(shí)驗(yàn)室委托測序中心進(jìn)行RNA-seq測序后得到的

12、最原始的數(shù)據(jù)為了便于測序數(shù)據(jù)的發(fā)布和共享，高通量測序數(shù)據(jù)以FASTQ格式來記錄所測的堿基讀段和質(zhì)量分?jǐn)?shù)如圖1所示，F(xiàn)ASTQ格式以測序讀段為單位存儲，每條讀段占4行，其中第1行和第3行由文件識別標(biāo)志和讀段名(ID)組成(第1行以“”開頭而第3行以“ ”開頭；第3行中ID可以省略，但“ ”不能省略)，第2行為堿基序列，第4行為對應(yīng)的測序質(zhì)量分?jǐn)?shù)關(guān)于FASTQ格式更多地介紹可參考文獻(xiàn)21為方便保存和共享各實(shí)驗(yàn)室產(chǎn)生/doc/be6b0724c324a6024ee96d8b的高通量測序數(shù)據(jù)，NCBI、EBI、DDBJ等數(shù)據(jù)中心建立了大容量的數(shù)據(jù)庫SRA(SequenceReadArchive,/Tr

13、aces/sra)來存放共享的測序數(shù)據(jù)22-23836生物化學(xué)與生物物理進(jìn)展Prog.Biochem.Biophys.2010;37(8)每4行標(biāo)識為一個測序讀段讀段識別碼堿基序列此處圖片未下載成功此處圖片未下載成功此處圖片未下載成功 讀段識別碼測序質(zhì)量分?jǐn)?shù)鄄seq數(shù)據(jù)的基本處理seq的基本應(yīng)用是測量一個樣本中的基因表達(dá)或轉(zhuǎn)錄組有實(shí)驗(yàn)表明，新一代高通量測序技術(shù)重復(fù)數(shù)據(jù)之間的相關(guān)度較高(R2抑0.96)1-2，因此，如果對同一樣本在多個通道上進(jìn)行了RNA測序的技術(shù)重復(fù)，我們建議可以把幾個通道的數(shù)據(jù)進(jìn)行合并，這樣等效地增加了測序深度本節(jié)討論單個樣本RNA測序數(shù)據(jù)的基本處理流程，如圖2a所示圖1讀段

14、FASTQ數(shù)據(jù)格式示例Fig.1FASTQformatexamples(a)(b)圖2所示的流程,虛線箭頭表示可選輸入.鄄seq數(shù)據(jù)處理和分析流程圖鄄seqdataprocessingandanalysis(a)RNA-seq數(shù)據(jù)的基本處理,其方法介紹見正文第3節(jié).(b)兩類樣本RNA-seq數(shù)據(jù)比較分析的框架,對應(yīng)于正文的第4節(jié).(b)中虛線框內(nèi)為(a)援1讀段定位獲得RNA-seq的原始數(shù)據(jù)后，首先需要將所有測序讀段通過序列映射(mapping)定位到參考基因組上，這是所有后續(xù)處理和分析的基礎(chǔ)在讀段定位之前，有時還需要根據(jù)測序數(shù)據(jù)情況對其做某些基本的預(yù)處理例如，過濾掉測序質(zhì)量較差的讀/do

15、c/be6b0724c324a6024ee96d8b段、對miRNA測序讀段數(shù)據(jù)去除接頭序列等高通量測序的海量數(shù)據(jù)對計算機(jī)算法的運(yùn)行時間提出了很高的要求針對諸如Illumina/Solexa等測序平臺得到的讀段一般較短、且插入刪除錯誤較;37(8)王曦等：新一代高通量RNA測序數(shù)據(jù)的處理與分析837少等特點(diǎn)，人們開發(fā)了一些短序列定位算法這些算法主要采用空位種子索引法(spaced-seedindexing)或Burrows-Wheeler轉(zhuǎn)換(Burrows-WheelerTransform，BWT)技術(shù)來實(shí)現(xiàn)24空位種子索引法首先將讀段切分，并選取其中一段或幾段作為種子建立搜索索引，再通過查

16、找索引、延展匹配來實(shí)現(xiàn)讀段定位，通過輪換種子考慮允許出現(xiàn)錯配(mismatch)的各種可能的位置組合BWT方法通過B-W轉(zhuǎn)換25將基因組序列按一定規(guī)則壓縮并建立索引，再通過查找和回溯來定位讀段，在查找時可通過堿基替代來實(shí)現(xiàn)允許的錯配表1列出了目前可免費(fèi)下載使用的部分短序列定位軟件其中采用空位種子片段索引法的代表是Maq26，而采用Burrows-Wheeler轉(zhuǎn)換的代表是Bowtie27總的來說，采用BWT的定位算法在時間效率上要優(yōu)于空位種子片段索引法24,28隨著讀長的增加，允許讀段序列中存在插入刪除(indel)的定位變得可行而重要由于以上兩類方法對序列中插入刪除的處理較為困難，近來人們開

17、發(fā)了一些基于改進(jìn)的Smith-Waterman動態(tài)規(guī)劃算法29的序列比對工具，如BFAST30、SHRiMP31、Mosaik(/marthlab/Mosaik)等，但算法速度較慢，大多需采用計算機(jī)并行編程技術(shù)來解決運(yùn)行時間的問題表1Table1名稱MAQ26Bowtie27BWA32ZOOM33ELANDSOAP234RazerS35NovoalignSHRiMP31BFAST30Mosaik)否是是否否/doc/be6b0724c324a6024ee96d8b否否是否質(zhì)量2)是是是否否否否是是適用于Illumina/Solexa測序平臺的讀段定位軟件Mappers/alignersforI

18、llumina/Solexasequencingdata主要采用技術(shù)空位種子BWTBWT空位種子空位種子BWTq-grams過濾Needleman-Wunsch算法空位種子q-grams過濾Smith-Waterman算法:/shrimp:/index.php/BFAST/marthlab/Mosaik網(wǎng)址:/:/index.shtml/bwa.shtml:/products/zoom:/software/genome_analyzer_software.ilmn:/www.seqan.de/projects/razers.html是是是是Waterman算法并行編程Smith-Waterma

19、n算法并行編程):是否能以SAM格式輸出;2)質(zhì)量:是否提供讀段定位質(zhì)量信息;BWT:Burrows-Wheeler轉(zhuǎn)換.在RNA測序數(shù)據(jù)的基因組定位中，一個特殊的問題是跨越兩個外顯子接合區(qū)的讀段(junctionreads)定位在真/doc/be6b0724c324a6024ee96d8b核生物中，成熟的mRNA是經(jīng)過由mRNA前體中的外顯子經(jīng)過剪接形成的如果一個讀段跨越了兩個外顯子，那么就無法將這個讀段完整地定位到基因組序列上而同時，這種跨兩個外顯子的讀段在分析轉(zhuǎn)錄本的剪接形式和研究選擇性剪接中有重要的作用為了解決這一問題，人們采取兩種典型的策略來進(jìn)行接合區(qū)讀段的定位：一是根據(jù)已知的基因外

20、顯子注釋，構(gòu)建所有可能的外顯子接合區(qū)序列，與基因組序列一并作為定位的參考基因組；二是不依賴基因注釋，而是先利用能完整定位到基因組的讀段得到粗略的外顯子區(qū)域，并結(jié)合剪接位點(diǎn)序列構(gòu)建出可能的剪接位點(diǎn)，然后將不能完整定位的讀段分段定位到兩個外顯子可能的結(jié)合區(qū)域Illumina/Solexa平臺提供的RNA-seq軟件分析包GApipeline采用了第一種策略采用第二種策略的軟件有Tophat36和G-Mo.R-Se37等，最新的Tophat軟件增加了利用已知外顯子邊界注釋信息的選項(xiàng)838生物化學(xué)與生物物理進(jìn)展不論是哪種測序平臺，測序中都不可避免地存在一定的錯誤，基因組中又存在單核苷酸多態(tài)性等引起的序

21、列變化，所以在讀段定位時通常允許一定數(shù)量的錯配，可以根據(jù)不同應(yīng)用調(diào)節(jié)允許錯配的程度另一方面，由于基因組中重復(fù)序列和高相似度序列的影響，某些讀段會出現(xiàn)定位到基因組多個位置的情況這些因素影響了各個讀段到基因組的定位質(zhì)量，在一些新的讀段定位算法中，同時給出每個讀段與基因組匹配質(zhì)量通常在后續(xù)處理前，人們將多定位的讀段都過濾掉，也有人嘗試用適當(dāng)?shù)牟呗园讯喽ㄎ蛔x段“分配”到其中某些位置上2,38.讀段定位到基因組后通常采用SAM(SequenceAlignment/Map)格式或其二進(jìn)制版本BAM格式39來存儲二進(jìn)制版本可大大節(jié)省存儲空間，但不能直接用普通文本編輯工具顯示關(guān)于SAM格式的詳細(xì)介紹，可查閱(

22、/SAM1.pdf)援2基因表達(dá)水平估計在深度測序技術(shù)出現(xiàn)之前，高通量測量不同基因表達(dá)水平的主要手段是基因芯片，在此基礎(chǔ)上可以對不同組織或者不同發(fā)育階段的基因表達(dá)差異和模式進(jìn)行分析mRN/doc/be6b0724c324a6024ee96d8bA-seq數(shù)據(jù)最基本的應(yīng)用也是檢測基因的表達(dá)水平，與基因芯片數(shù)據(jù)相比，RNA測序得到的是數(shù)字化的表達(dá)信號，具有靈敏度高、分辨率高、無飽和區(qū)等優(yōu)勢40-42測序數(shù)據(jù)是對提取出的RNA轉(zhuǎn)錄本中隨機(jī)進(jìn)行的短片段測序，如果一個轉(zhuǎn)錄本的豐度高，則測序后定位到其對應(yīng)的基因組區(qū)域的讀段也就多，可以通過對定位到基因外顯子區(qū)的讀段計數(shù)來估計基因表達(dá)水平很顯然，讀段計數(shù)除了

23、與基因真實(shí)表達(dá)水平成正比，還與基因長度成正比，同時也與測序深度即測序?qū)嶒?yàn)中得到的總讀段數(shù)正相關(guān)為了保持對不同基因和不同實(shí)驗(yàn)間估計的基因表達(dá)值的可比性，人們提出了RPM和RPKM的概念2RPM(readspermillionreads)即每百萬讀段中來自于某基因的讀段數(shù)，考慮了測序深度對讀段計數(shù)的影響RPKM(readsperkilobasespermillionreads)是每百萬讀段中來自于某基因每千堿基長度的讀段數(shù)，公式表示為：=基因區(qū)讀段計數(shù)伊伊109不僅對測序深度作了歸一化，而且對基.Biochem.Biophys.2010;37(8)因長度也作了歸一化，使得不同長度的基因在不同測序深

24、度下得到的基因表達(dá)水平估計值具有了可比性，是目前最常用的基因表達(dá)估計方法軟件rSeq43、DEGseq軟件包44和Cufflinks45等都提供了用上述方法進(jìn)行基因表達(dá)水平計算的功能根據(jù)RNA-seq文庫制備標(biāo)準(zhǔn)，在不考慮基因結(jié)構(gòu)的理想情況下，讀段會均勻地分布在基因上而實(shí)際上，通過對實(shí)際數(shù)據(jù)的可視化分析很容易發(fā)現(xiàn)，讀段在基因上的分布有著自身的一些模式，呈現(xiàn)出不均勻性(圖3)這一問題已經(jīng)引起很多學(xué)者的關(guān)注46-48造成讀段分布出現(xiàn)偏好的原因可能有多個方面：在制備cDNA文庫時，反轉(zhuǎn)錄所采用的隨機(jī)引物對RNA序列具有一定的偏好性，使得cDNA片段不能夠完全均勻地取自各轉(zhuǎn)錄本；在PCR擴(kuò)增中，擴(kuò)增效

25、率與序列的GC含量等特征相關(guān)，可導(dǎo)致GC含量高的cDNA片段在文庫中拷貝數(shù)增加超過其他片段；舍棄多定位的讀段也可能導(dǎo)致讀段的非均勻分布；等等如果能對讀段分/doc/be6b0724c324a6024ee96d8b布的不均勻性進(jìn)行建模并加以校正，可以提高RNA-seq推斷基因表達(dá)量的準(zhǔn)確度但根據(jù)對實(shí)際數(shù)據(jù)的觀察，對于較長轉(zhuǎn)錄本，讀段非均勻分布帶來的誤差很大程度上可相互抵消，用RPKM來估計基因的表達(dá)水平可以得到比較滿意的結(jié)果援3選擇性剪接事件識別和剪接異構(gòu)體表達(dá)水平推斷在真核生物中，選擇性剪接現(xiàn)象普遍存在基因轉(zhuǎn)錄形成的mRNA前體(pre-mRNA)在剪接過程中因去掉不同的內(nèi)含子區(qū)域或保留不同的

26、外顯子區(qū)域，可形成不同的剪接異構(gòu)體根據(jù)RNA-seq原理，只要測序深度足夠深，就能檢測到所有轉(zhuǎn)錄本的全部序列，包括來自剪接接合區(qū)的序列利用考慮到接合區(qū)的讀段定位方法，就有可能系統(tǒng)地研究某一組織或某一條件下的基因選擇性剪接事件前面已經(jīng)提到，Tophat等軟件定位剪接接合區(qū)讀段的策略能標(biāo)定出剪接事件中的兩個剪接位點(diǎn)：供體位點(diǎn)和受體位點(diǎn)通過比較供體位點(diǎn)和受體位點(diǎn)的組合，就能識別選擇性剪接事件4,49圖3中包含了選擇性剪接識別的一個例子進(jìn)一步，通過對供體和受體位點(diǎn)的讀段計數(shù)，結(jié)合外顯子其他區(qū)域的讀段數(shù)據(jù)，還能定量地計算選擇性剪接事件之間的比例50-51對于每一個剪接異構(gòu)體，RNA-seq數(shù)據(jù)能在一;3

27、7(8)王曦等：新一代高通量RNA測序數(shù)據(jù)的處理與分析839定程度上推斷其表達(dá)水平比如，可以根據(jù)已知外顯子組成和各外顯子長度對剪接異構(gòu)體建立數(shù)學(xué)模型，在測序讀段轉(zhuǎn)錄本上均勻分布的假設(shè)下，利用各外顯子上的讀段數(shù)和接合區(qū)讀段數(shù)求解異構(gòu)體的表達(dá)值Jiang等43的方法及軟件IsoInfer52和cufflinks45都采用了這種思路來實(shí)現(xiàn)剪接異構(gòu)體的表達(dá)推斷需要指出的是，某些形式的剪接異構(gòu)體表達(dá)水平在這種方法框架下不可辨識533援4新基因的檢測在對RNA-seq數(shù)據(jù)的分析中，人們發(fā)現(xiàn)，往往不是所有讀段都能定位到已有注釋的基因區(qū)，說明除了轉(zhuǎn)錄噪聲或測序錯誤等的影響外，可能還存在尚未被注釋的基因這里，我

28、們把這種尚未注釋的基因稱為新基因，包括新的蛋白質(zhì)編碼基因和非編碼/doc/be6b0724c324a6024ee96d8bRNA基因能檢測新基因，尤其是低表達(dá)基因是RNA-seq技術(shù)優(yōu)于基因芯片的特點(diǎn)之一，因?yàn)樗恍枰靡阎蜃⑨寔碓O(shè)計檢測探針seq技術(shù)靈敏度高，但樣品污染、測序錯誤等仍可能帶來背景噪聲從基因組未注釋區(qū)域的RNA測序讀段信號中檢測新基因是典型的信號檢測問題如何控制新基因識別的誤發(fā)現(xiàn)率(FDR)是檢測方法的關(guān)鍵Useq軟件包54將ChIP-seq數(shù)據(jù)分析的方法移植到RNA-seq數(shù)據(jù)上，用滑窗的方法來識別測序讀段定位富集的區(qū)域，給出反映滑窗所在區(qū)域讀段富集顯著程度的P值(P-

29、value)及新基因誤發(fā)現(xiàn)率，通過設(shè)定P值或誤發(fā)現(xiàn)率的閾值，可篩選出讀段富集的區(qū)域，再將相鄰區(qū)域合并或根據(jù)剪接接合區(qū)讀段將相應(yīng)區(qū)域連接，完成新基因的檢測援5讀段的可視化及注釋對于復(fù)雜的組學(xué)數(shù)據(jù)，能盡可能方便地直接觀察數(shù)據(jù)對于數(shù)據(jù)的分析和解釋都非常重要，對新一代測序數(shù)據(jù)的可視化和交互展示是一個非常重要但容易被人忽視的問題不深入考查數(shù)據(jù)的細(xì)節(jié)，而是滿足于對數(shù)據(jù)的統(tǒng)計分析，是高通量數(shù)據(jù)應(yīng)用中經(jīng)常容易陷入的誤區(qū)，方便有效的可視化工具能夠幫助避免這樣的誤區(qū)表2列出了部分適用于RNA-seq數(shù)據(jù)的全基因組瀏覽器，其中比較具有代表性的有UCSCGenomeBrowser、CisGenomeBrowser和I

30、GV(IntegrativeGenomicsViewer)等這些瀏覽器具有如下特點(diǎn)：a能在不同尺度下顯現(xiàn)單個或多個讀段在基因組上的位置，包括來源于剪接接合區(qū)的讀段；b能在不同尺度下顯示不同區(qū)域的讀段豐度，以反映不同區(qū)域的轉(zhuǎn)錄水平或測序效率；c能顯示基因及其剪接異構(gòu)體的注釋信息；d能顯示其他注釋信息，例如物種間基因組序列保守性、序列GC含量等；e能直接或間接支持SAM/BAM讀段定位數(shù)據(jù)存儲格式.UCSCGenomeBrowser55屬于基于網(wǎng)絡(luò)模式的全基因組瀏覽器，所有數(shù)據(jù)都需要上載到遠(yuǎn)程服務(wù)器，經(jīng)過處理后將圖形返回客戶端顯示圖3中的例子就是從UCSCGenomeBrowser的顯示截取的.C

31、isGenomeBrowser56是典型的本地版基因組瀏覽器，所有讀段數(shù)據(jù)、注釋信息都存于本地文件，因此不/doc/be6b0724c324a6024ee96d8b需要網(wǎng)絡(luò)連接，方便內(nèi)部考查數(shù)據(jù)用IGV(/igv)可以說是以上兩種模式的融合，既可以從遠(yuǎn)程服務(wù)器端下載各種注釋信息，又可以從本地加載注釋信息表2名稱IGV適用于mRNA鄄seq數(shù)據(jù)的全基因組瀏覽器/viewersapplicabletomRNA鄄seqdataviewing支持的數(shù)據(jù)格式網(wǎng)址:/igv/GFF3,BED,SAM/BAM,WIG,55BED,bigBed,BEDGRAPH,GFF,GTF,bigWig,BAM,5657

32、,BED,refFlat,FA,Wig,BED,GFF,FASTA,ELAND,GFF,BED,MAQ,SAMBAM,BED,GFF2,GFF3,FASTA,VCFExpressiondata,Annotationtracks:/jiangh/browser/sj/mochiview-start:/www.bioinformatics.bbsrc.ac.uk/projects/seqmonk/p/gambit-viewer:/packages/release/bioc/html/GenomeGraphs.html以上列出的全基因組瀏覽器均可/doc/be6b0724c324a6024ee96d

33、8b在Windows、Linux和蘋果公司的MacOS等計算機(jī)平臺下運(yùn)行.840生物化學(xué)與生物物理進(jìn)展Prog.Biochem.Biophys.2010;37(8)剪接接合區(qū)樣本A測序標(biāo)簽分布該基因在A、B樣本中差異表達(dá)內(nèi)含子區(qū)的測序標(biāo)簽接測序標(biāo)布剪接接合區(qū)樣本B測序標(biāo)簽分布基因注釋圖3鄄seq數(shù)據(jù)可視化示例鄄seqdatavisualization圖示區(qū)域?yàn)槿祟惢駽BX7.圖中紅色表示樣本A的數(shù)據(jù)，藍(lán)色表示樣本B.各軌道(track)依次為:基因組坐標(biāo)、樣本A的剪接接合區(qū)、樣本A的讀段分布、樣本B的剪接接合區(qū)、樣本B的讀段分布、UCSC基因注釋.圖中還標(biāo)識了：因受體位點(diǎn)不同而形成的選擇性剪

34、接;基因的5憶端出現(xiàn)讀段的非均勻分布;在兩個樣本中,差異表達(dá)基因的讀段信號強(qiáng)度不同;在基因標(biāo)注的內(nèi)含子(intron)區(qū)域存在少量不連續(xù)的讀段.除對讀段的可視化外，用描述統(tǒng)計學(xué)方法對實(shí)驗(yàn)數(shù)據(jù)進(jìn)行分類別統(tǒng)計也十分重要例如，統(tǒng)計讀段在各個染色體上的分布情況和在注釋的外顯子、內(nèi)含子、剪接接合區(qū)、基因間區(qū)的分布情況等目前，已經(jīng)有一些用于測序數(shù)據(jù)注釋的生物信息學(xué)軟件，比如SAMtools39、BEDtools58等，但由于測序技術(shù)發(fā)展迅速，用戶需求因人而異，用戶經(jīng)常還需要根據(jù)需求編寫一定的程序或腳本完成或完善注釋分析的任務(wù)對于熟悉圖形用戶界面的研究人員，還可以利用UCSCTableBrowse56和Ga

35、laxy59-60來配合完成注釋分析由于UCSCTableBrowser集成了大量基因組尺度上的注釋信息，而Galaxy又為用戶提供了書寫簡單、接口明晰和直觀的數(shù)據(jù)處理流程，這種方法十分方便有效，也為很多學(xué)者在展示研究成/doc/be6b0724c324a6024ee96d8b果時所采用以上描述了對于基因組已知的物種進(jìn)行RNA-seq數(shù)據(jù)處理基本流程，圖2a給出了其主要步驟若研究對象尚未完成基因組測序，則需要采用讀段的從頭拼裝(denovoassembly)6,61-62來代替讀段定位，后續(xù)流程也須做相應(yīng)的調(diào)整若RNA-seq實(shí)驗(yàn)在文庫制備時保留了RNA的方向信息，則應(yīng)分別研究來自正鏈和反鏈的

36、轉(zhuǎn)錄產(chǎn)物，并通過與基因注釋比較來檢測反義轉(zhuǎn)錄本63最近，RNA-MATE64軟件在其分析流程中加入如此的處理策略此外，通過分析定位到外顯子接合區(qū)的讀段，還可以獲取轉(zhuǎn)錄本結(jié)構(gòu)，這為研究基因的剪接調(diào)控機(jī)理提供了重要信息5而利用RNA-seq數(shù)據(jù)提供的序列信息，通過與DNA序列的細(xì)致比較可分析轉(zhuǎn)錄組的序列差異(如SNP等)65，從而研究等位基因的表達(dá)模式66-67及RNA編輯68等最后需要指出，由于miRNA在序列和結(jié)構(gòu)上具有一定的特點(diǎn)，miRNA-seq數(shù)據(jù)的基本處理流程也與本節(jié)所述有所不同，感興趣的讀者可參考軟件工具miRDeep69提供的處理策略多類樣本mRNA鄄seq數(shù)據(jù)間的比較分析很多RN

37、A-seq實(shí)驗(yàn)的目的是為了比較兩種或多種樣本中基因表達(dá)或整個轉(zhuǎn)錄組的差異，如比較癌癥組織和正常組織的轉(zhuǎn)錄組差異等這些差異既包括通常意義下的差異表達(dá)基因，也主要包括選擇性剪接模式的差異、剪接異構(gòu)體表達(dá)的差異、非編碼轉(zhuǎn)錄本的差異等這些差異一般可以用一些統(tǒng)計假設(shè)檢驗(yàn)方法檢測，但這種檢驗(yàn)有時會受到測序深度、基因長度等因素的影響70-71，需要對結(jié)果進(jìn)行仔細(xì)分析，消除可能的混雜因素，必要時可以用讀段的絕對表達(dá)值倍數(shù)變化(fold-change)來作為補(bǔ)充圖2b給出了兩類樣本數(shù)據(jù)分析的框架;37(8)王曦等：新一代高通量RNA測序數(shù)據(jù)的處理與分析841援1差異表達(dá)基因的識別雖然新一代測序相對第一代測序的單

38、位成本大大降低，但是，利用RNA測序進(jìn)行基因表達(dá)研究的成本仍很高，因此，很多實(shí)驗(yàn)室沒有條件進(jìn)行樣本重復(fù)如果兩類樣本均沒有生物重復(fù)，例如只/doc/be6b0724c324a6024ee96d8b對兩個細(xì)胞系各進(jìn)行一次mRNA樣本測序，則可以用隨機(jī)采樣模型通過假設(shè)檢驗(yàn)來分析差異表達(dá)對于某個基因，如果一個讀段來自于這個基因，我們稱事件A發(fā)生對于一次RNA-seq實(shí)驗(yàn)，事件A發(fā)生的概率可以用這個基因上的讀段數(shù)n除以所有基因上的讀段總數(shù)N來估計，即RPM事件A發(fā)生的概率反應(yīng)了這個基因的表達(dá)水平如果要判斷(a)某個基因在兩個樣本中的表達(dá)水平是否一致，就可以通過檢驗(yàn)事件A在兩種條件下發(fā)生的概率是否一致來實(shí)

39、現(xiàn)，采用似然比檢驗(yàn)1、Fisher精確檢驗(yàn)72以及基于MA圖的統(tǒng)計檢驗(yàn)方法44等同樣，也可用RPKM作為統(tǒng)計量來進(jìn)行假設(shè)檢驗(yàn)分析，由于是比較同一個基因在兩個樣本間的差異，基因長度的影響被抵消，用RPKM和用RPM得到的結(jié)果相似對無生物重復(fù)的RNA-seq數(shù)據(jù)進(jìn)行差異表達(dá)基因分析，已經(jīng)有幾個公開發(fā)表的軟件，包括DEGseq44、Useq54、Cufflinks45中的Cuffdiff模塊等圖4展示了我們開發(fā)的DEGseq軟件提供的多種差異表達(dá)基因識別方法的應(yīng)用例子(b)-A(d)A-0A8log2(readcountsforeachgene)inB(c)-AvsBAvsB圖4用DEGseq軟件包

40、識別差異表達(dá)基因的結(jié)果(a)各基因在樣本A和樣本B中表達(dá)水平的散點(diǎn)圖.(b),(c),(d)圖中紅點(diǎn)表示分別用FET、LRT和MARS方法得到的差異表達(dá)基因.FET:FishersExactTest,Fisher精確檢驗(yàn).LRT:LikelihoodRatioTest,似然比檢驗(yàn).MARS:MA-plot-basedmethodwit/doc/be6b0724c324a6024ee96d8bhRandomSamplingmodel,基于MA圖的隨機(jī)采樣模型.如果每一類樣本都包含了若干生物重復(fù)，如病人和正常人對照研究，則可以沿用基因芯片數(shù)據(jù)分析中的很多方法比如，可以用t檢驗(yàn)結(jié)合倍數(shù)變化的方法來分

41、析差異表達(dá)如果兩類樣本具有配對的信息，也可以通過整合每對樣本分析結(jié)果來實(shí)現(xiàn)其步驟為，先在每對樣本中識別出差異表達(dá)的基因，再尋找這若干組差異表達(dá)基因之間的相同者，或用投票的方法來為基因的差異程度打分針對某些RNA-seq數(shù)據(jù)生物樣本量小，R軟件包DEGseq44和edgeR73等還專門提供了基于改進(jìn)模型的統(tǒng)計方法此外，一類將分類器與特征選擇包裹在一起的方法也同樣適用于此類問題(見4.3)842生物化學(xué)與生物物理進(jìn)展4援2差異表達(dá)剪接異構(gòu)體的識別差異表達(dá)剪接異構(gòu)體的識別方法與差異表達(dá)基因的識別相似如果把剪接異構(gòu)體看成是獨(dú)立的基因，那么前面討論的用于識別差異表達(dá)基因的方法對剪接異構(gòu)體完全適用但是，注

42、意到來自于同一個基因的剪接異構(gòu)體并不獨(dú)立，某些假設(shè)檢驗(yàn)的基本條件并不滿足，得到的結(jié)果就不一定正確此外，由于現(xiàn)在剪接異構(gòu)體表達(dá)推斷的方法還不夠成熟，加之在基因結(jié)構(gòu)不可辨識的剪接異構(gòu)體上作表達(dá)推斷會出現(xiàn)病態(tài)結(jié)果，差異表達(dá)剪接異構(gòu)體的識別問題還處于探索的階段目前，在剪接異構(gòu)體表達(dá)水平可辨識且讀段覆蓋度較高的基因上，BASIS74方法通過貝葉斯模型來推斷差異表達(dá)的剪接異構(gòu)體換一個角度，剪接異構(gòu)體由選擇性剪接造成，如果剪接異構(gòu)體的表達(dá)有差異，那么導(dǎo)致這些異構(gòu)體的選擇性剪接事件及異構(gòu)體特異的外顯子的表達(dá)也會有差異因此，對差異表達(dá)剪接異構(gòu)體的識別可以轉(zhuǎn)變?yōu)榉治鲞x擇性剪接事件和外顯子表達(dá)的差異75外顯子表達(dá)差

43、異的分析可以完全利用基因表達(dá)差異的分析方法而剪接接合區(qū)也可以看成是一個較短的“外顯子”(長度一般與測序長度相當(dāng))不過，由于外顯子長度較基因的長度短，對應(yīng)的讀段數(shù)量較少，差異識別的敏感度會有所下降Solas方法75就是根據(jù)類似的原理，采用統(tǒng)計學(xué)假設(shè)檢驗(yàn)的方法來識別差異表達(dá)的剪接異構(gòu)體4援3對樣本的分類分析:/doc/be6b0724c324a6024ee96d8b通過統(tǒng)計方法識別出來的差異表達(dá)基因及剪接異構(gòu)體能否有效地區(qū)別兩類樣本，可以通過分類分析進(jìn)一步證實(shí)如果把每個基因(或剪接異構(gòu)體)的表達(dá)值作為特征，則差異表達(dá)基因(或剪接異構(gòu)體)的選取也就是特征篩選的過程把前面用統(tǒng)計方法等檢測出來的差異表達(dá)

44、基因(或剪接異構(gòu)體)用于分類分析，常被稱為過濾法另一類基于分類器的包裹法，例如R-SVM76、SVM-RFE77等，可以根據(jù)每個特征在分類器中所占的權(quán)重來篩選特征，因此也可以用于差異表達(dá)基因(或剪接異構(gòu)體)的識別分類的性能可以用交叉驗(yàn)證(cross-validation，CV)方法來評估需要特別注意的是，交叉驗(yàn)證應(yīng)該包括對特征選擇步驟的交叉驗(yàn)證，防止發(fā)生信息泄露而導(dǎo)致評估結(jié)果過于樂觀具體做法是：將樣本按一定的策略分成兩份，一份(通常是樣本數(shù)多的一份)用于特征選取和分類器訓(xùn)練，而用余下的樣本進(jìn)行分類器性能的估計；重復(fù)以上步驟多次，就得.Biochem.Biophys.2010;37(8)到交叉驗(yàn)

45、證錯誤率必要時還可以用隨機(jī)置換檢驗(yàn)(permutationtest)來推斷所得錯誤率的統(tǒng)計顯著性76當(dāng)樣本數(shù)較小時，可以采用留一法交叉驗(yàn)證(leave-one-outcross-validation，LOOCV)4援4其他高層分析方法檢測差異表達(dá)的基因或差異表達(dá)異構(gòu)體是人們認(rèn)識所研究的生物問題機(jī)理的第一步，接下來需要從功能上研究這些差異轉(zhuǎn)錄現(xiàn)象的分子機(jī)理這與在基因芯片應(yīng)用中所面臨的是同樣的生物學(xué)問題，對芯片數(shù)據(jù)分析結(jié)果的后續(xù)處理方法，都可以借鑒到測序數(shù)據(jù)上來如何進(jìn)一步地從機(jī)理來解釋結(jié)果，還需結(jié)合已知生物學(xué)知識進(jìn)行后續(xù)分析人們對基因芯片得到的基因表達(dá)數(shù)據(jù)進(jìn)行分析的很多方法都可以用到RNA-seq

46、數(shù)據(jù)上來，比如利用機(jī)器學(xué)習(xí)方法進(jìn)行分類和特征選擇，對差異表達(dá)的基因進(jìn)行GO(geneontology)78類別富集分析、信號通路富集分析等，一些常用的分析工具包括GoMiner79、DAVID80和VisANT81等需要說明的是，在各種以差異表達(dá)基因?yàn)榛A(chǔ)的分析中，由于基因表達(dá)水平都是通過讀段計數(shù)來估計的，表達(dá)水平較高或/doc/be6b0724c324a6024ee96d8b轉(zhuǎn)錄本較長的基因擁有更多的讀段，更容易被多數(shù)統(tǒng)計方法識別為差異表達(dá)基因70這種偏好可能對后續(xù)分析帶來影響以GO類別富集分析為例，這種偏好將導(dǎo)致長基因占主導(dǎo)的功能類別更有可能被識別為富集的功能這將對生物機(jī)理的研究帶來誤導(dǎo)最

47、近，Young等71發(fā)展了一種GOseq方法，針對這一偏好對GO類別富集分析做了改進(jìn)RNA鄄seq數(shù)據(jù)處理中的生物信息學(xué)挑戰(zhàn)高通量測序技術(shù)的發(fā)展十分迅速，這要求相應(yīng)的數(shù)據(jù)處理與分析方法快速跟進(jìn)正是這些方法，架起了高通量實(shí)驗(yàn)數(shù)據(jù)與科學(xué)問題之間的橋梁這種橋梁作用正日趨重要，也為生物信息學(xué)帶來了挑戰(zhàn)7,71這里，我們重點(diǎn)討論兩方面的挑戰(zhàn)：a如何實(shí)現(xiàn)剪接接合區(qū)讀段的準(zhǔn)確定位？b在數(shù)據(jù)處理各階段中，如何對RNA-seq數(shù)據(jù)的系統(tǒng)誤差和固有偏好建?；蜓a(bǔ)償，以消除它們可能帶來的錯誤推斷及結(jié)論？援1剪接接合區(qū)讀段的定位測序技術(shù)的一個發(fā)展趨勢是測序長度不斷增加隨著讀長的增加，RNA-seq中來自剪接接合區(qū)的讀段

48、會越來越多我們粗略估算，按照人類基因組refSeq基因注釋，一般情況下，如果測序讀長;37(8)王曦等：新一代高通量RNA測序數(shù)據(jù)的處理與分析843為50個堿基，則約有10%的讀段來自剪接接合區(qū)而當(dāng)測序長度達(dá)到100個堿基時，這個比例將達(dá)到25%左右對這些剪接接合區(qū)讀段的分析，將使我們能夠更準(zhǔn)確地檢測剪接事件和推斷剪接異構(gòu)體的表達(dá)水平，大大推進(jìn)人們對選擇性剪接的研究在RNA-seq出現(xiàn)的早期，人們沒有意識到剪接接合區(qū)讀段的重要性因?yàn)楫?dāng)時的讀長只有2030個堿基，來自剪接接合區(qū)的讀段所占比例甚小當(dāng)時讀段定位的通常做法是，先將讀段與全基因組序列做映射定位，再考慮不能定位的讀段是否來自于剪接接合區(qū)2

49、這種做法雖然在一定程度上保證了讀段定位的比率，但由于基因組中重復(fù)序列和相似序列的存在，部分接合區(qū)讀段有可能在容許錯配的情況下被定位到基因組上其他位置，從而失去了定位到正確的剪接接合區(qū)的機(jī)會在讀段定位時，如果要同時/doc/be6b0724c324a6024ee96d8b考慮基因組序列和剪接接合區(qū)序列，就要利用已知的剪接事件注釋，這是目前軟件通用的方法然而，包括人類在內(nèi)的各物種的基因注釋信息都還有待完善，也沒有較完整的剪接組(splicome)數(shù)據(jù)庫，能夠不依賴注釋信息和對剪接機(jī)理的現(xiàn)有認(rèn)識，高效、準(zhǔn)確地定位所有已知和未知的接合區(qū)讀段，仍然是對讀段映射定位算法的一個挑戰(zhàn)援2系統(tǒng)噪聲和偏好的分析雖

50、然深度測序技術(shù)的準(zhǔn)確性較以前的技術(shù)有了很大提高，但仍然存在錯誤和噪聲比如從圖3中可以看到，內(nèi)含子區(qū)內(nèi)有一些不連續(xù)的讀段，很可能由系統(tǒng)噪聲造成，如樣品污染、測序錯誤和不恰當(dāng)?shù)淖x段定位策略等從圖3還能看出，外顯子區(qū)域內(nèi)的讀段信號分布也很不均勻有文獻(xiàn)報道，序列組成尤其是GC含量46、RNA二級結(jié)構(gòu)2等也有可能是導(dǎo)致讀段不均勻分布的原因這些噪聲和分布偏好將影響新基因的識別和對剪接異構(gòu)體形式和表達(dá)水平的推斷合理地建模RNA-seq數(shù)據(jù)中的系統(tǒng)噪聲和偏好是解決上述問題最有效的辦法基本的思路可以是：首先根據(jù)實(shí)驗(yàn)原理尋找可能產(chǎn)生系統(tǒng)噪聲或偏差的因素，并盡可能將這些因素轉(zhuǎn)化成可量化的特征，如序列特征、二級結(jié)構(gòu)等

51、；然后，將用實(shí)驗(yàn)數(shù)據(jù)對這些特征做統(tǒng)計分析，構(gòu)造和訓(xùn)練模型，用模型來對數(shù)據(jù)進(jìn)行校正需要注意的是，某些偏好是由當(dāng)前的測序技術(shù)和分析方法共同造成的，難以完全消除71在這種情況下，后續(xù)處理和解釋時需要充分意識到這種偏好可能對生物學(xué)結(jié)論帶來的影響，必要時通過補(bǔ)充其他實(shí)驗(yàn)來驗(yàn)證和修正通過高通量測序得到的生物學(xué)結(jié)論總結(jié)與展望本文以Illumina/Solexa測序平臺為例，嘗試對新一代測序技術(shù)的RNA-seq數(shù)據(jù)處理和分析方法做了較為全面的梳理，并對各個環(huán)節(jié)上可用的軟件進(jìn)行了匯總高通量測序是正在飛速發(fā)展的技術(shù)，相應(yīng)的生物信息學(xué)方法也在快速發(fā)展，這里討論的是RNA-seq中一些代表性的方法和問題，希望能對正在

52、或即將采用RNA-seq實(shí)驗(yàn)進(jìn)行科學(xué)研究的學(xué)者和進(jìn)行RNA測序數(shù)據(jù)處理的同行提供參考.測序和基因芯片有很多共同的應(yīng)用領(lǐng)域，盡管相對還不是很成熟，RNA-seq技術(shù)在很多方面已經(jīng)表現(xiàn)/doc/be6b0724c324a6024ee96d8b出了優(yōu)勢，有人甚至預(yù)言基因芯片時代即將結(jié)束36但也有報道認(rèn)為，RNA-seq數(shù)據(jù)在基因表達(dá)水平的估計上和基因芯片相比沒有明顯的優(yōu)勢60，加上測序的成本目前還遠(yuǎn)高于芯片實(shí)驗(yàn)的成本，所以更多人認(rèn)為測序和基因芯片將長期共存，以各自不同的特點(diǎn)在現(xiàn)代組學(xué)研究中發(fā)揮作用新一代高通量測序技術(shù)的應(yīng)用面非常廣82，RNA-seq只是其中一個方面，除此之外，基因組的從頭測序和重測

53、序83-84、染色質(zhì)免疫沉淀測序(ChIP-seq)85-86、甲基化測序(Methyl-seq)87-88等技術(shù)都同樣有著廣泛的應(yīng)用尤其是，用ChIP-seq研究蛋白質(zhì)與DNA的相互作用，能夠得到高分辨率的轉(zhuǎn)錄因子結(jié)合數(shù)據(jù)和組蛋白修飾等表觀遺傳學(xué)數(shù)據(jù)發(fā)展有效的生物信息學(xué)方法，將ChIP-seq數(shù)據(jù)與RNA-seq得到的轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行綜合分析，將大大推進(jìn)人們對復(fù)雜的基因轉(zhuǎn)錄調(diào)控系統(tǒng)的認(rèn)識致謝感謝本實(shí)驗(yàn)室劉霖曦、謝芃、孟璐等同學(xué)對本工作有意義的討論，感謝斯坦福大學(xué)WingHWong教授、HuiJiang博士和JunLi同學(xué)等的討論和幫助參考文獻(xiàn)1MarioniJC,MasonCE,ManeSM,

54、etal.RNA-seq:anassessment.GenomeRes,2008,18(9):1509-15172MortazaviA,WilliamsBA,McCueK,etal.MappingandSeq.NatMethods,2008,5(7):621-628844生物化學(xué)與生物物理進(jìn)展3NagalakshmiU,WangZ,WaernK,etal.Thetranscriptional:/doc/be6b0724c324a6024ee96d8bing.Science,2008,320(5881):1344-13494SultanM,SchulzMH,RichardH,etal.Aglob

55、alviewofgene.Science,2008,321(5891):956-9605WangET,SandbergR,LuoS,etal.Alternativeisoformregulation.Nature,2008,456(7221):470-4766BirzeleF,SchaubJ,RustW,etal.Intotheunknown:expression.NucleicAcidsRes,2010,doi:10.1093/nar/gkq1167SangerF,NicklenS,CoulsonAR.DNAsequencingwithchain-terminatinginhibitors.

56、ProcNatlAcadSciUSA,1977,74(12):5463-54678MarguliesM,EgholmM,AltmanWE,etal.Genomesequencingdensitypicolitrereactors.Nature,2005,437(7057):376-3809ShendureJ,JiH.Next-generationDNAsequencing.NatBiotechnol,26(10):1135-114510RuparelH,BiL,LiZ,etal.Designandsynthesisofa3憶-O-ally/doc/be6b0724c324a6024ee96d8

57、bl.ProcNatlAcadSciUSA,2005,102(17):5932-593711SeoTS,BaiX,KimDH,etal.Four-colorDNAsequencingby.ProcNatlAcadSciUSA,2005,102(17):5926-593112JuJ,KimDH,BiL,etal.Four-colorDNAsequencingbysynthesis.ProcNatlAcadSciUSA,2006,103(52):19635-1964013FedurcoM,RomieuA,WilliamsS,etal.BTA,anovelreagentforphaseamplifi

58、edDNAcolonies.NucleicAcidsRes,2006,34(3):e2214ShendureJA,PorrecaGJ,ChurchGM.OverviewofDNA/AusubelFM,BrentR,KingstonRE,A:JohnWileyandSons,Inc.,2008:Unit7.115MardisER.Next-generationDNAsequencingmethods.AnnuRev,2008,9:387-40216Full/doc/be6b0724c324a6024ee96d8berCW,MiddendorfLR,BennerSA,etal.Thechallen

59、gesof.NatBiotechnol,2009,27(11):1013-102317CroucherNJ,FookesMC,PerkinsTT,etal.Asimplemethodfor.NucleicAcidsRes,2009,37(22):e14818ParkhomchukD,BorodinaT,AmstislavskiyV,etal.TranscriptomespecificsequencingofcomplementaryDNA.NucleicAcidsRes,2009,37(18):e12319PerkinsTT,KingsleyRA,FookesMC,etal.Astrand-s

60、pecificSeqanalysisofthetranscriptomeofthetyphoidbacillusSalmonellatyphi.PLoSGenet,2009,5(7):e100056920MamanovaL,AndrewsRM,JamesKD,etal.FRT-seq:free,strand-specifictranscriptomesequencing.Nat.Biochem.Biophys.;37(8),2010,7(2):130-13221CockPJ,FieldsCJ,GotoN,etal.TheSangerFASTQfileformat,andtheSolexa/Il