RNA-seq數(shù)據(jù)分析
我們對(duì)工具和技術(shù)的評(píng)估落腳點(diǎn)在它們鑒定出的差異基因的準(zhǔn)確性端辱。為了完成這個(gè)評(píng)估,至少需要四個(gè)不同的分析階段(圖2;表2)鹅经。第一階段把測(cè)序平臺(tái)生成的原始測(cè)序數(shù)據(jù)比對(duì)到轉(zhuǎn)錄組锅很。第二階段量化與每個(gè)基因或轉(zhuǎn)錄本來源的reads數(shù)量窿冯,構(gòu)建表達(dá)矩陣赵抢。該過程可能包括1個(gè)或多個(gè)子過程如比對(duì)剧蹂,組裝和定量,或者它也可以一個(gè)從讀取計(jì)數(shù)生成表達(dá)矩陣烦却。通常有一個(gè)第三階段宠叼,包括過濾低表達(dá)的基因和至關(guān)重要的移除樣品間技術(shù)差異的標(biāo)準(zhǔn)化過程。DGE的最后階段是構(gòu)建樣本分組和其它協(xié)變量的統(tǒng)計(jì)模型其爵,計(jì)算差異表達(dá)置信度冒冬。
第1階段-測(cè)序reads的比對(duì)和組裝
測(cè)序完成后,分析的起點(diǎn)是包含測(cè)序堿基的FASTQ文件摩渺。最常見的第一步是將測(cè)序reads比對(duì)到已知的轉(zhuǎn)錄組(或注釋的基因組)窄驹,將每個(gè)測(cè)序reads轉(zhuǎn)換為一個(gè)或多個(gè)基因組坐標(biāo)。傳統(tǒng)上证逻,該過程是通過幾個(gè)不同的比對(duì)工具(如TopHat,STAR或HISAT)完成的,其都依賴參考基因組的存在囚企。由于測(cè)序的cDNA來自RNA丈咐,可能跨越外顯子邊界,因此與參考基因組(包含內(nèi)含子和外顯子)比對(duì)時(shí)需要進(jìn)行剪接比對(duì)龙宏,即允許reads中出現(xiàn)大片段gap棵逊。
如果沒有可用的包含已知外顯子邊界的高質(zhì)量基因組注釋,或者如果希望將reads與轉(zhuǎn)錄本(而不是基因)相關(guān)聯(lián)银酗,則需要在比對(duì)后執(zhí)行轉(zhuǎn)錄組組裝步驟辆影。諸如StringTie和SOAPdenovo-Trans之類的組裝工具使用比對(duì)reads的gap來推測(cè)外顯子邊界和可能的剪接位點(diǎn)。轉(zhuǎn)錄本重頭組裝特別適用于參考基因組注釋缺失或不完整的物種黍特,或者對(duì)異常轉(zhuǎn)錄本感興趣(例如在腫瘤組織中)的研究蛙讥。轉(zhuǎn)錄組組裝方法受益于雙端測(cè)序和/或更長(zhǎng)的reads的使用,增加跨越splice junctions的可能性灭衷。但是次慢,通常不需要從RNA-seq數(shù)據(jù)中從頭做轉(zhuǎn)錄組組裝來確定DGE (生信寶典注:無參分析組裝是必須的)。
最近翔曲,涌現(xiàn)了一些計(jì)算效率高的“alignment free”工具迫像,例如Sailfish,Kallisto和Salmon瞳遍,它們將測(cè)序reads直接與轉(zhuǎn)錄本關(guān)聯(lián)闻妓,而無需單獨(dú)的定量步驟。這些工具在定量高豐度(以及長(zhǎng)度更長(zhǎng))的轉(zhuǎn)錄本方面表現(xiàn)出很好的性能掠械。但是由缆,它們?cè)诙康拓S度或短轉(zhuǎn)錄本方面不夠準(zhǔn)確。(39個(gè)工具份蝴,120種組合深度評(píng)估 (轉(zhuǎn)錄組分析工具哪家強(qiáng)))
第2階段-定量轉(zhuǎn)錄本豐度
將reads比對(duì)到基因組或轉(zhuǎn)錄組后犁功,下一步就是將它們分配給基因或轉(zhuǎn)錄本,獲得表達(dá)矩陣婚夫。不同的比較研究表明浸卦,定量過程中采用的方法對(duì)最終結(jié)果的影響最大,甚至比比對(duì)工具影響更大案糙。
常用的定量工具包括RSEM限嫌,CuffLinks,MMSeq和HTSeq时捌,以及上述的無比對(duì)直接定量工具怒医。基于reads計(jì)數(shù)的工具(例如HTSeq或featureCounts)通常會(huì)丟棄許多比對(duì)的序列奢讨,包括那些具有多個(gè)匹配位置或比對(duì)到多個(gè)表達(dá)特征的reads稚叹。這可以在隨后的分析中消除同源和重疊的轉(zhuǎn)錄本。RSEM使用期望最大化模型來分配模糊的reads,而無參考的比對(duì)方法(例如Kallisto)則將這些reads用于后續(xù)的定量扒袖,這可能會(huì)導(dǎo)致結(jié)果偏差塞茅。轉(zhuǎn)錄本豐度估計(jì)可以轉(zhuǎn)換成等效的read計(jì)數(shù),能完成這一轉(zhuǎn)換的部分工具依賴tximport包季率。量化步驟結(jié)束后會(huì)得到一個(gè)合并的表達(dá)矩陣野瘦,每個(gè)表達(dá)特征(基因或轉(zhuǎn)錄本)各占一行,每個(gè)樣品各占一列飒泻,中間的值是實(shí)際讀數(shù) (reads count)或估計(jì)的表達(dá)豐度鞭光。
階段3-過濾和標(biāo)準(zhǔn)化
通常,基因或轉(zhuǎn)錄本的reads count需要進(jìn)行過濾和標(biāo)準(zhǔn)化泞遗,以移除測(cè)序深度惰许、表達(dá)模式和技術(shù)偏差的影響。過濾去除在所有樣本中都低豐度表達(dá)的基因是很直接的方式刹孔,標(biāo)準(zhǔn)化表達(dá)矩陣的方法要復(fù)雜一些啡省。簡(jiǎn)單的轉(zhuǎn)換可以校正豐度,降低GC含量和測(cè)序深度的影響
階段4-差異表達(dá)分析
獲得表達(dá)矩陣后髓霞,就可以構(gòu)建統(tǒng)計(jì)模型評(píng)估哪些轉(zhuǎn)錄本發(fā)生了顯著的表達(dá)改變卦睹。有幾個(gè)常用工具可以完成此任務(wù);一些基于基因水平的表達(dá)計(jì)數(shù)方库,其它的基于轉(zhuǎn)錄本水平的表達(dá)計(jì)數(shù)结序。基因水平的工具通常依賴于比對(duì)的reads計(jì)數(shù)纵潦,并使用廣義線性模型來進(jìn)行復(fù)雜實(shí)驗(yàn)設(shè)計(jì)的評(píng)估徐鹤。這些工具包括EdgeR,DESeq2和limma + voom等工具邀层,這些工具計(jì)算效率高并且彼此之間結(jié)果穩(wěn)定性好返敬。評(píng)估差異isoforms表達(dá)的工具,例如CuffDiff寥院,MMSEQ和Ballgown劲赠,往往需要更多的計(jì)算資源,并且結(jié)果的變化也更大秸谢。但是凛澎,在差異表達(dá)工具應(yīng)用之前的操作(即關(guān)于比對(duì)、定量估蹄、過濾和標(biāo)準(zhǔn)化)對(duì)最終結(jié)果的影響更大塑煎。
參考基因組和基因注釋文件獲取
通常測(cè)序生成的reads要與參考基因組或參考轉(zhuǎn)錄組進(jìn)行比對(duì),或Pseudo-alignment臭蚁。所以首先需要獲取參考基因組和參考轉(zhuǎn)錄組信息最铁。
Ensembl http://www.ensembl.org/info/data/ftp/index.html 是常用的信息齊全的參考基因組和GTF文件下載網(wǎng)站讯赏。
Ensembl提供的參考基因組有2種組裝形式和3種重復(fù)序列處理方式, 分別是primary, toplevel和unmasked (dna)、soft-masked (dna_sm)和masked (dna_rm)炭晒。一般選擇dna.primary或dna_sm.primary野建。
為什么選擇Primary
Primary assembly contains all toplevel sequence regions excluding haplotypes and patches. This file is best used for performing sequence similarity searcheswhere patch and haplotype sequences would confuse analysis.
為什么不選擇masked
Masked基因組是指所有重復(fù)區(qū)和低復(fù)雜區(qū)被N代替的基因組序列榴啸,比對(duì)時(shí)就不會(huì)有reads比對(duì)到這些區(qū)域。一般不推薦用masked的基因組艳汽,因?yàn)樗斐闪诵畔⒌膩G失嗤无,由此帶來的一個(gè)問題是uniquely比對(duì)到masked基因組上的reads實(shí)際上可能不是unique的震束。而且masked基因組還會(huì)帶來比對(duì)錯(cuò)誤,使得在允許錯(cuò)配的情況下当犯,本來來自重復(fù)區(qū)的reads比對(duì)到基因組的其它位置垢村。 另外檢測(cè)重復(fù)區(qū)和低復(fù)雜區(qū)的軟件不可能是完美的,這就造成遮蓋住的重復(fù)序列和低復(fù)雜區(qū)并不一定是100%準(zhǔn)確和敏感的嚎卫。
soft-masked基因組是指把所有重復(fù)區(qū)和低復(fù)雜區(qū)的序列用小寫字母標(biāo)出的基因組嘉栓,由于主要的比對(duì)軟件,比如BWA拓诸、bowtie2等都忽略這些soft-mask侵佃,直接把小寫字母當(dāng)做大寫字母比對(duì),所以使用soft-masked基因組的比對(duì)效果和使用unmasked基因組的比對(duì)效果是相同的奠支。
基因注釋GTF文件在分析轉(zhuǎn)錄組數(shù)據(jù)時(shí)會(huì)用到馋辈,也從這獲取,GTF文件的解釋見文件格式部分倍谜。ENSEMBL的基因注釋文件與GeneCode(http://www.gencodegenes.org/)V26版本一致
下載基因功能和結(jié)構(gòu)注釋信息
ENSEMBL數(shù)據(jù)庫(kù)的BioMart http://www.ensembl.org/biomart/martview工具為下載基因的功能信息迈螟、序列信息、結(jié)構(gòu)信息尔崔、ID的轉(zhuǎn)換等提供了很大的便利答毫。
注意在BioMart的Attribute選項(xiàng)里如果選擇了蛋白相關(guān)的選項(xiàng),得到的結(jié)果中只有蛋白編碼基因的信息季春。如果要下載所有基因信息洗搂,請(qǐng)不要選擇蛋白相關(guān)的選項(xiàng)。
GTF/GFF文件格式解讀和轉(zhuǎn)換
測(cè)序原始數(shù)據(jù)下載
生物或醫(yī)學(xué)中涉及高通量測(cè)序的論文鹤盒,一般會(huì)將原始測(cè)序數(shù)據(jù)上傳到公開的數(shù)據(jù)庫(kù)蚕脏,上傳方式見測(cè)序文章數(shù)據(jù)上傳找哪里;并在文章末尾標(biāo)明數(shù)據(jù)存儲(chǔ)位置和登錄號(hào)
NCBI的SRA (Sequence Read Archive) 數(shù)據(jù)庫(kù)(http://www.ncbi.nlm.nih.gov/sra/) 是最常用的存儲(chǔ)測(cè)序數(shù)據(jù)的數(shù)據(jù)庫(kù)侦锯。目前SRA數(shù)據(jù)的組織方式分為下面4個(gè)層次:
Studies—研究課題驼鞭;
Experiments—實(shí)驗(yàn)設(shè)計(jì);
Runs—測(cè)序結(jié)果集尺碰;
Samples—樣品信息挣棕。
使用NCBI提供的SRA-toolkit中的工具fastq-dump直接下載SRR文件译隘,并轉(zhuǎn)換為FASTQ格式,--split-3參數(shù)表示如果是雙端測(cè)序就自動(dòng)拆分洛心,如果是單端不受影響固耘。--gzip轉(zhuǎn)換fastq為壓縮文件,節(jié)省空間词身。