转录组数据分析my RNA-seqRNASeq 数据分析

简单小需求:如何将FPKM转换成TPM?

2019-02-21  本文已影响10人  刘小泽

刘小泽写于19.2.21

或许这还是个比较常用的需求

我们做转录组分析,得到的结果可能是FPKM。但是FPKM确实存在不准确性,推荐使用TPM。

read count和FPKM结果都可以转成TPM,但是因为FPKM跟TPM的计算都考虑了基因长度,所以从FPKM转TPM最方便快捷。

假设原来的表达矩阵fpkm_expr.txt中行为基因,列为样本,中间数值是FPKM计算得到的值

先读取自己的表达矩阵

expMatrix<-read.table("fpkm_expr.txt",header = T,row.names = 1)

其实早已经有人帮我们整理好了

https://haroldpimentel.wordpress.com/2014/05/08/what-the-fpkm-a-review-rna-seq-expression-units/

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)
}

如果要计算TPM值,只需要用一下apply函数

tpms <- apply(expMatrix,2,fpkmToTpm)
tpms[1:3,]

最后可以根据TPM的特征进行检查,看每列加和是否一致

colSums(tpms)
上一篇下一篇

猜你喜欢

热点阅读