在高通量測序數(shù)據(jù)的分析中，僅僅靠raw read counts描述基因的表達(dá)量是遠(yuǎn)遠(yuǎn)不夠的牺荠。受限于測序過程中的技術(shù)因素影響翁巍，read counts對于基因表達(dá)量的反映存在一定偏好（bias）。因此休雌，Mortazavi等人提出了RPKM/FPKM的方法對read counts進(jìn)行normalization以使基因表達(dá)量的比較可以在不同文庫間進(jìn)行灶壶。隨后，Mortazavi等人更是提出了考慮轉(zhuǎn)錄本長度分布情況的TPM方法杈曲。本文將會簡要說明為什么我們要對read counts進(jìn)行normalization驰凛，以及RPKM，F(xiàn)PKM担扑，TPM是什么洒嗤，并通過一個簡單地例子闡述為什么TPM才是被更多人認(rèn)同的方法。
為什么我們要進(jìn)行Normalization
我們之所以說測序得到的read counts并不是其mRNA豐度的忠實反映魁亦，是因為read counts的數(shù)量會受到多種因素的影響渔隶，例如：
-
測序深度：某些低表達(dá)量的基因只有在較高的測序深度時才能檢測到。一般而言洁奈，隨著測序深度的增加间唉，基因種類以及可變剪接體的數(shù)目也會增加。同時利术，測序深度高的樣本read counts也會較高呈野。
Effect of sequencing depth
在上圖中，樣本A中的基因表達(dá)量是樣本B的兩倍印叁，但這是由于測序深度引起的結(jié)果被冒，而非真實存在的差異军掂。
-
基因長度：由于高通量測序是將cDNA碎片化后再進(jìn)行測序的，因此越長的基因產(chǎn)生的碎片也會更多昨悼，在測序中也會更加容易被檢測到蝗锥。所以對基因長度的校正也是十分有必要的。
Effect of gene length
在這個圖里率触，基因X和Y的真實表達(dá)量是一致的终议，但是基因X的reads會比基因Y要多，這是由于基因X的基因長度較長所致的葱蝗。
除了上述兩個主要因素外穴张，還會有其他因素對read counts的檢測有所影響，例如轉(zhuǎn)錄組的組成两曼，GC含量皂甘，random hexamers引起的測序偏好等等。由于上述因素的存在悼凑，導(dǎo)致在不同樣本間使用read counts 進(jìn)行比較是不太現(xiàn)實的偿枕，人們便提出了許多對read counts進(jìn)行Normalization的方法。本文在此簡單地介紹使用最為廣泛卻最受質(zhì)疑的RPKM/FPKM佛析，以及受到更多人認(rèn)可的TPM益老。
RPKM/FPKM
RPKM: Reads per kilo base per million mapped reads
FPKM: Fragments per kilo base per million mapped reads
從RPKM的公式中我們可以看到，RPKM對基因長度（gene length）以及測序深度（mapped reads from library）都進(jìn)行了校正寸莫。而FPKM只是RPKM的雙端測序（pair-end）版本捺萌。
TPM
TPM: Transcript per million
TPM的計算公式如下所示：
同RPKM一樣，TPM對基因的長度進(jìn)行了校正膘茎，計算比對到基因上的reads/基因長度得到長度校正的表達(dá)量 reads per kilobase (RPK)桃纯。再以文庫中RPK之和作為Scale Factor求出TPM。
相比于RPKM使用read counts之和來作為文庫校正因子披坏，TPM使用RPK之和作為文庫校正因子的好處是考慮了不同樣本間的基因長度的分布态坦。因為RPK是一個對基因長度進(jìn)行校正后的表達(dá)量單位，所以RPK之和也不會再帶入基因長度的bias棒拂。因此伞梯，如果需要比較的樣本之間轉(zhuǎn)錄本分布不一致時（例如不同物種RNA-seq的比較），使用TPM是一個較佳的Normalization方案帚屉。
一個例子
以下有ABCD四個基因谜诫，并同時進(jìn)行了三次重復(fù)的測序（rep1，rep2攻旦，rep3）喻旷。首先，可以看出由于基因長度的關(guān)系牢屋，B基因的read counts都是較其余的基因較高的（排除D而言）且预。其次槽袄，可以看出rep3的測序深度較高，得到的read counts也較多锋谐，同時還檢測到其余兩個重復(fù)中沒有檢測到的D基因（低豐度基因）遍尺。不管是測序深度還是基因長度的嚴(yán)重地干擾了我們對不同樣本的基因表達(dá)量比較。因此怀估，我們采用RPKM和TPM的方法進(jìn)行Normalization后狮鸭，再來比較其中的差異合搅。
為了方便起見多搀，以下的函數(shù)省略了kilo base，million base的轉(zhuǎn)換
RPKM <- function(x, gene.length){
lib.size <- sum(x)
rpm <- x/lib.size
rpkm <- rpm/gene.length
return(rpkm)
}
TPM <- function(x, gene.length){
rpk <- x/gene.length
scale.factor <- sum(rpk)
tpm <- rpk/scale.factor
return(tpm)
}
expr <- data.frame(rep1=c(10,20,5,0), rep2=c(12,25,8,0), rep3=c(30,60,15,1),
row.names = c("A", "B", "C", "D"))
gl <- c(2,4,1,10)
rpkm1 <- apply(expr, 2, FUN = RPKM, gene.length = gl)
tpm1 <- apply(expr, 2, FUN = TPM, gene.length = gl)
# 對結(jié)果都乘1000灾部，看起來直觀點
rbind(rpkm1*1000, colSums(rpkm1)*1000)
## rep1 rep2 rep3
## A 142.8571 133.3333 141.5094340
## B 142.8571 138.8889 141.5094340
## C 142.8571 177.7778 141.5094340
## D 0.0000 0.0000 0.9433962
## 428.5714 450.0000 425.4716981
rbind(tpm1*1000, colSums(tpm1)*1000)
## rep1 rep2 rep3
## A 333.3333 296.2963 332.594235
## B 333.3333 308.6420 332.594235
## C 333.3333 395.0617 332.594235
## D 0.0000 0.0000 2.217295
## 1000.0000 1000.0000 1000.000000
通過以上的例子康铭，我們可以發(fā)現(xiàn)RPKM和TPM最重要且直接的區(qū)別就是，TPM校正后的表達(dá)量之和在不同樣本中是相等的赌髓，而在RPKM中并不是从藤。這就導(dǎo)致當(dāng)我們使用RPKM在樣本間進(jìn)行比較的時候就不能得到可靠的結(jié)果。
雖然以上的例子較為理想化锁蠕，但在進(jìn)行實驗時夷野，我們認(rèn)同的一個前提是不同的技術(shù)重復(fù)（例如一個細(xì)胞重復(fù)測三次）所得到的基因表達(dá)量應(yīng)當(dāng)是一致的。在測序的過程中或許會引入一定的誤差荣倾，但這些誤差是可以通過Normalization的方法在一定程度上校正回來的悯搔。簡而言之，RPKM和TPM這類方法就是為了使不同樣本間的總體表達(dá)量趨于一致舌仍，讓不同樣本間的基因表達(dá)量有可比較性妒貌，而TPM在這方面表現(xiàn)得較為出色，故受到學(xué)界更多的青睞铸豁。
- 常用的Normalization 方法總結(jié)
| Normalization method | 描述 | 考慮因素 | 使用場景 |
|---|---|---|---|
| CPM (counts per million) | 使用read counts的總和校正counts | 測序深度 | 同一樣本組重復(fù)樣本之間的基因counts比較灌曙；不適用于樣品內(nèi)的比較或差異分析 |
| TPM (transcripts per kilobase million) | 每百萬mapped reads中每kb轉(zhuǎn)錄本上的reads數(shù) | 測序深度和基因長度 | 樣本內(nèi)或同一樣本組樣本之間的基因counts比較； 不適用于差異分析 |
| RPKM/FPKM (reads/fragments per kilobase of exon per million reads/fragments mapped) | 如TPM | 測序深度和基因長度 | 樣本內(nèi)基因間的counts比較; 不適用于樣本間比較和差異分析 |
| DESeq2’s median of ratios [1] | counts除以樣本特異的校正因子节芥，該因子由基因計數(shù)相對于每個基因的幾何平均值的中位數(shù)比率確定 | 測序深度及RNA組成 | 樣本之間的基因counts比較以及差異分析; 不適用于樣本內(nèi)比較 |
| EdgeR’s trimmed mean of M values (TMM) [2] | 使用樣本之間的加權(quán)截尾的對數(shù)表達(dá)量比值的均值進(jìn)行TMM校正 | 測序深度在刺，RNA組成以及基因長度 | 樣品之間和樣品內(nèi)部的基因counts比較，適用于差異分析 |
這邊補充一點是在差異表達(dá)分析中头镊，對基因表達(dá)量的Normalization沒有這么嚴(yán)謹(jǐn)蚣驼，只是要求對文庫大小進(jìn)行校正即可，例如
edgeR包中的cpm(counts per million)只是簡單地將read counts/library size拧晕。對這方面感興趣的同學(xué)可以看我介紹
edgeR的那篇文章：簡單的轉(zhuǎn)錄組基因差異表達(dá)分析 -- edgeR其中有更詳細(xì)的介紹隙姿。
參考文章:
- https://www.rna-seqblog.com/rpkm-fpkm-and-tpm-clearly-explained/
- https://www.biostars.org/p/273537/
- http://www.reibang.com/p/1940c5954c81
- Mortazavi A, Williams BA, McCue K, Schaeffer L, Wold B. Mapping and quantifying mammalian transcriptomes by RNA-Seq. Nature methods. 2008 Jul 1;5(7):621-8.
- https://hbctraining.github.io/DGE_workshop/lessons/02_DGE_count_normalization.html
完。
