TPM格式搜索错了，好尴尬

Posted on by

网上的答案经常不靠谱

通常情况下我会使用 featureCounts 得到表达矩阵是 raw counts, 但总是有人需要我转换成各种形式，比如 RPKM, FPKM and TPM

概念见： https://statquest.org/2015/07/09/rpkm-fpkm-and-tpm-clearly-explained/

想偷一下懒，就搜索到了 这个答案： http://ny-shao.name/2016/11/18/a-short-script-to-calculate-rpkm-and-tpm-from-featurecounts-output.html

## functions for rpkm and tpm
## from https://gist.github.com/slowkow/c6ab0348747f86e2748b#file-counts_to_tpm-r-L44
## from https://www.biostars.org/p/171766/
rpkm <- function(counts, lengths) {
  rate <- counts / lengths
  rate / sum(counts) * 1e9
}

tpm <- function(counts, lengths) {
  rate <- counts / lengths
  rate / sum(rate) * 1e6
}

朋友提醒我，其实错了，因为很明显 colSums(exprSet_tpm) 并不是一百万。

其实我老早就写过TPM公式，就是RPKM的百分比的百万倍扩大值，所以还是自己动手重新写了代码。

rpkm <- function(counts, lengths) {
  rate <- counts / lengths
  rate / sum(counts) * 1e9
}

exprSet_rpkm=rpkm(exprSet,len) 
exprSet_tpm=1e6*exprSet_rpkm/colSums(exprSet_rpkm)

不过，比较奇怪的是这个时候 colSums(exprSet_tpm) 是接近一百万，而不是精确的一百万，我还没有想清楚具体缘由，是不是R的计算小数点问题。

有关于的讨论： http://blog.nextgenetics.net/?e=51#body-anchor

Comments are closed.