美文网首页my RNA-seq
网上的答案经常不靠谱

网上的答案经常不靠谱

作者: 因地制宜的生信达人 | 来源:发表于2018-07-01 20:51 被阅读137次

    通常情况下我会使用 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

    相关文章

      网友评论

        本文标题:网上的答案经常不靠谱

        本文链接:https://www.haomeiwen.com/subject/bmyhuftx.html