Counts FPKM RPKM TPM CPM 的转化

最近有粉丝自告奋勇希望可以把他自己在简书等平台的生物信息学笔记分享在我们生信技能树公众号，在专业的舞台上跟大家切磋！

非常欢迎啊！

下面是简书的昵称:  嘿嘿嘿嘿哈  ，的分享

一、什么是RPKM、 FPKM、TPM、CPM

RPKM, FPKM and TPM, clearly explained - StatQuest!!!
对基因counts进行校正定量一般有RPKM、FPKM 、TPM和CPM这几种方法，StatQuest网站中对RPKM, FPKM 和 TPM作了通俗简要的说明，现将其核心要点整理如下：

RPKM

RPKM (Reads Per Kilobase Million, or Reads Per Kilobase of transcript per Million reads mapped)
每千个碱基的转录每百万映射读取的reads

计算RPKM的方法：

1. 计算样本中的总reads，并将该数字除以1,000,000——这是我们的“每百万”缩放因子 ( “per million” scaling factor ) 。

2. 将reads数除以“每百万”缩放因子。消除测序深度影响，得到每百万reads（RPM, reads per million ）

3. 将RPM值除以基因长度（以千碱基为单位），消除基因长度影响，得到RPKM。

FPKM

FPKM (Fragments Per Kilobase Million, or Fragments Per Kilobase of transcript per Million reads mapped) 每千个碱基的转录每百万映射读取的fragments

FPKM与RPKM非常相似。RPKM是针对单端测序的RNA-seq而言的，其中每个reads对应于一个已测序的单个片段。FPKM用于双端测序的RNA-seq。使用双端测序RNA-seq，两个reads可以对应一个片段（Fragment）。RPKM和FPKM之间的唯一区别是FPKM考虑到两次reads可以映射到一个片段（因此它不会对该片段进行两次计数）。
即 单端测序：reads=fragments，双端测序：2 * reads≈fragments
而经过上游处理，双端测序两个reads可以对应一个片段的过程已经完成，最后得到的counts就已经相当于是片段fragments了，因此下游分析由counts计算RPKM、 FPKM这两者的公式完全一致。

TPM

TPM (Transcripts Per Million, or Transcripts Per kilobase of exon model per Million mapped reads)
每千个碱基的转录每百万映射读取的Transcripts

计算TPM的方法：

1. 将读数计数除以每个基因的长度（以千碱基为单位），得到每千碱基reads（RPK, reads per kilobase）。

2. 计算样本中所有RPK值，然后将其除以1,000,000，得到“每百万”缩放因子 ( “per million” scaling factor )。

3. 将RPK值除以“每百万”缩放因子，得到TPM。

因此在计算TPM时，与RPKM、 FPKM相比，唯一的区别是TPM先对基因长度进行标准化，然后对测序深度进行标准化。但是，这种差异的影响非常深远。因为使用TPM时，每个样本中所有TPM的总和是相同的，这样可以更轻松地比较每个样本中映射到基因的读数的比例。相反，使用RPKM和FPKM，每个样本中的标准化读数之和可能会有所不同，这使得直接比较样本变得更加困难。

CPM

CPM（Counts Per Million, or Counts of exon model Per Million mapped reads） 每百万映射读取的counts

除了RPKM、 FPKM、TPM这几种方法，CPM也是较为常见的一种基因定量方式。原始的表达量除以该样本表达量的总和，再乘以一百万，即可得到CPM值。CPM值只对测序深度进行了标准化，一般利用edgeR包的cpm()函数即可对基因counts进行简单校正 。
edgeR::cpm(counts)

二、由Counts计算FPKM/RPKM和TPM

有许多文章已经给出了这几种计数方式的计算和转化关系，如What the FPKM? A review of RNA-Seq expression units | The farrago (wordpress.com)。

countToTpm <- function(counts, effLen) 
{ 
  rate <- log(counts) - log(effLen) denom <- log(sum(exp(rate))) exp(rate - denom + log(1e6)) 
} 
countToFpkm <- function(counts, effLen) { 
  N <- sum(counts) exp( log(counts) + log(1e9) - log(effLen) - log(N) ) 
} 
fpkmToTpm <- function(fpkm) { 
  exp(log(fpkm) - log(sum(fpkm)) + log(1e6)) 
} 
countToEffCounts <- function(counts, len, effLen) { 
  counts * (len / effLen) 
} 
################################################################################ 
# An example 
################################################################################ 
cnts <- c(4250, 3300, 200, 1750, 50, 0) 
lens <- c(900, 1020, 2000, 770, 3000, 1777) 
countDf <- data.frame(count = cnts, length = lens)
# assume a mean(FLD) = 203.7 
countDf$effLength <- countDf$length - 203.7 + 1 
countDf$tpm <- with(countDf, countToTpm(count, effLength)) 
countDf$fpkm <- with(countDf, countToFpkm(count, effLength)) 
with(countDf, all.equal(tpm, fpkmToTpm(fpkm))) 
countDf$effCounts <- with(countDf, countToEffCounts(count, length, effLength))

为了更清楚直观地明白FPKM/RPKM和TPM的计算过程，根据StatQuest定义给出的计算方法，自己重新写了由counts转化为FPKM/RPKM和TPM的计算公式代码，经过验证其结果与上面公式所得结果一致，那今后就用自己写的这些代码吧，具体代码如下：

#TPM (Transcripts Per Kilobase Million)  每千个碱基的转录每百万映射读取的Transcripts 
counts2TPM <- function(count=count, efflength=efflen){   
  RPK <- count/(efflength/1000)   #每千碱基reads (reads per kilobase) 长度标准化   
  PMSC_rpk <- sum(RPK)/1e6        #RPK的每百万缩放因子 (“per million” scaling factor ) 深度标准化   
  RPK/PMSC_rpk                       
}                    
#FPKM/RPKM (Fragments/Reads Per Kilobase Million )  每千个碱基的转录每百万映射读取的Fragments/reads 
#RPKM与FPKM分别针对单端与双端测序而言，计算公式是一样的 
counts2FPKM <- function(count=count, efflength=efflen){    
  PMSC_counts <- sum(count)/1e6   #counts的每百万缩放因子 (“per million” scaling factor) 深度标准化   
  FPM <- count/PMSC_counts        #每百万reads/Fragments (Reads/Fragments Per Million) 长度标准化   
  FPM/(efflength/1000)                                       
} 
#FPKM与TPM的转化 
FPKM2TPM <- function(fpkm){   
  fpkm/sum(fpkm)*1e6   
} 
tpm <- as.data.frame(apply(counts,2,counts2TPM)) 
colSums(tpm) 

fpkm <- as.data.frame(apply(counts,2,counts2FPKM)) 
colSums(fpkm) 

tpm0 <- as.data.frame(apply(counts,2,FPKM2TPM)) 
colSums(tpm0)

要注意一点的是计算FPKM/RPKM和TPM时，基因长度一般指的都是基因的有效长度effective length，即该基因的外显子总长度或转录本总长度，以此为标准来消除测序造成的基因长度影响才更为准确。

可参见生信技能树文章：基因长度之多少 | 生信菜鸟团 (bio-info-trainee.com)

关于如何得到基因有效长度，请参阅文章：获取基因有效长度的N种方法。

参考资料：RPKM, FPKM and TPM, clearly explained - StatQuest!!!What the FPKM? A review of RNA-Seq expression units | The farrago (wordpress.com)RPKM、FPKM、TPM详解RNA-seq的counts值，RPM, RPKM, FPKM, TPM 的异同 - 云+社区 - 腾讯云 (tencent.com)

文末友情宣传

强烈建议你推荐给身边的博士后以及年轻生物学PI，多一点数据认知，让他们的科研上一个台阶：

• 数据挖掘（GEO,TCGA,单细胞）2022年6月场，快速了解一些生物信息学应用图表
• 生信入门课-2022年6月场，你的生物信息学第一课