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)而生了。

image

<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所示）。

image

<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)單，最直白的語(yǔ)言描述“geneA的FPKM是10”的測(cè)序意義。