https://zhuanlan.zhihu.com/p/50811365
Hello大家好!好久不見了吉嫩!
之前手頭上一直有很多事情,因此咱們的生物信息學(xué)100個(gè)基礎(chǔ)問題(BBQ100)也耽誤了一陣子唤殴,給大家鞠躬道歉斋攀,以后希望能夠保持一定的更新速度,早日填完我們這個(gè)立下的Flag弧哎!
根據(jù)之前的規(guī)劃雁比,我們將用接下來的幾期問題來探索一下RNA-Seq定量的問題,也就是要探索一下我們常說的RPKM撤嫩,F(xiàn)PKM偎捎,TPM,raw count 和RSEM,前面4個(gè)指標(biāo)都比較直觀茴她,方便理解寻拂,最后一個(gè)RSEM需要涉及到一些機(jī)器學(xué)習(xí)的知識(shí),我們盡量給大家把比較復(fù)雜的問題簡(jiǎn)單化丈牢,方便大家的入門祭钉。
1. RNA-Seq定量過程中的比較問題
我們?cè)?a target="_blank">BBQ-34的時(shí)候討論過RNA-Seq的方法論相關(guān)的問題,就是RNA-Seq的基本假設(shè)是什么己沛?簡(jiǎn)單來說就是 細(xì)胞/組織/個(gè)體 的兩種不同狀態(tài)進(jìn)行比較慌核,比較的目的就是尋找差異表達(dá)gene,然后從差異表達(dá)gene來推斷造成生理狀態(tài)不同的原因申尼。
而我們的RNA-Seq一般情況下是針對(duì)mRNA以及帶polyA的lncRNA進(jìn)行建庫(kù)測(cè)序分析的垮卓。那么理論上把測(cè)序的FASTQ文件mapping到參考基因組上,再結(jié)合參考基因組的GTF/GFF文件就可以找到全基因組的每一個(gè)gene上mapping到了多少個(gè)reads count师幕。
拿到了reads count以后粟按,我們就會(huì)嘗試著想要比較gene之間的表達(dá)量的關(guān)系,但是這時(shí)候往往會(huì)面臨兩個(gè)問題霹粥,舉個(gè)例子:
- 問題1: 比如我有g(shù)ene3灭将,有1000條測(cè)序reads,gene4有2000條測(cè)序reads蒙挑,那么我能否說gene4就一定比gene3的表達(dá)量高宗侦?(圖1 gene3 與 gene4)
- 問題2: 比如我有g(shù)ene1,有1000條測(cè)序reads忆蚀,我的另一個(gè)處理?xiàng)l件下gene2有2000條測(cè)序reads矾利,我能否就說geneA在處理?xiàng)l件下表達(dá)量降低了?(圖1 gene1與gene2)
在面臨這些比較問題的時(shí)候馋袜,我們就需要對(duì)mapping到gene的reads count進(jìn)行矯正男旗,至少根據(jù)問題1我們知道應(yīng)該在矯正的時(shí)候考慮過gene長(zhǎng)度的問題;根據(jù)問題2欣鳖,我們大概應(yīng)該能夠猜想到察皇,矯正的時(shí)候應(yīng)該需要考慮整體測(cè)序量的問題。到此泽台,RPKM和FPKM這兩個(gè)指標(biāo)就應(yīng)運(yùn)而生了什荣。
<figcaption style="margin-top: 0.66667em; padding: 0px 1em; font-size: 0.9em; line-height: 1.5; text-align: center; color: rgb(153, 153, 153);">圖1 ( Manuel Garber et al., Nature Methods, 2011 )</figcaption>
2. 什么是RPKM與FPKM?
RPKM = Reads Per Kilobase per Million mapped reads
假設(shè)回貼到geneA 的 reads count為 CountA怀酷,geneA的exon總長(zhǎng)度為L(zhǎng)en(A) Kbp稻爬,總的測(cè)序量為D兆reads,那么:
geneA RPKM = CountA / Len(A) / D * 10^9
那么什么是FPKM呢蜕依?先來看一下FPKM的定義:
FPKM = Fragments Per Kilobase per Million mapped reads
大家可以比較清楚看出來桅锄,RPKM中的R指的是Reads琉雳,F(xiàn)PKM中的F是指Fragments,Reads都比較好理解友瘤,就是我們的測(cè)序短的片段翠肘,那么fragment是什么呢?這是以為我們現(xiàn)在測(cè)序一般來說都是測(cè)雙端測(cè)序(paired-end sequencing)辫秧,那么在mapping回參考基因組的時(shí)候就會(huì)有兩條reads束倍,分別是read1和read2,分別來源于建庫(kù)打斷的5' 端和3'端茶没。那么這2條reads就可以在參考基因組上確定1個(gè)小的片段肌幽,這個(gè)片段就叫fragment(圖2所示)晚碾。
<figcaption style="margin-top: 0.66667em; padding: 0px 1em; font-size: 0.9em; line-height: 1.5; text-align: center; color: rgb(153, 153, 153);">圖2 (Frances S. Turner)</figcaption>
所以抓半,如果是現(xiàn)在最常用的雙端測(cè)序,1個(gè)gene的FPKM應(yīng)該等于RPKM / 2格嘁。
3. RPKM / FPKM有什么優(yōu)缺點(diǎn)笛求?
因?yàn)楝F(xiàn)在使用Illumina測(cè)序平臺(tái),絕大多數(shù)的測(cè)序都是使用雙端測(cè)序糕簿,那么基本上我們一般對(duì)gene進(jìn)行定量都是使用FPKM來進(jìn)行探入。FPKM的優(yōu)點(diǎn)大家都很了解了,能夠矯正掉gene長(zhǎng)度以及測(cè)序深度對(duì)gene表達(dá)定量的影響懂诗,那么FPKM的缺點(diǎn)大家是否熟悉呢蜂嗽?
一個(gè)比較容易被人提及的問題是對(duì)于不同批次測(cè)序的結(jié)果,所有g(shù)ene的FPKM的總和不是一個(gè)固定的值殃恒。比如WT 測(cè)的所有g(shù)ene的FPKM總和可能是10000植旧,treat組測(cè)到的FPKM總和可能是15000,這樣對(duì)于WT和treat組之間的差異表達(dá)gene的尋找就有可能出現(xiàn)問題离唐,這個(gè)時(shí)候就需要用到我們常用的另一種矯正方法TPM病附。
4. 提問環(huán)節(jié)
好了,相信通過今天的介紹亥鬓,大家能夠?qū)PKM與RPKM有一個(gè)比較清楚的認(rèn)識(shí)了完沪。我做一個(gè)簡(jiǎn)單的小提問:請(qǐng)用最簡(jiǎn)單,最直白的語言描述“geneA的FPKM是10”的測(cè)序意義嵌戈。