TCGA和GTEx的数据联合分析实战
0.缘起
很多文章中用到GEPIA这个网页工具来进行TCGA和GTEx表达量的联合比较,但在此之前我不知道要怎样在R语言中实现。这个GEPIA的文章里说:
The imbalance between the tumor and normal data can cause inefficiency in various differential analyses. Fortunately, the GTEx project produced RNA-Seq data for over 8000 normal samples, albeit from unrelated donors. Such data cannot be directly combined for integrated analysis due to many differences in aspects like data processing pipelines and gene models. To make data from different sources more compatible, the UCSC Xena project (http://xena.ucsc.edu/) has recomputed all expression raw data based on a standard pipeline to** minimize differences from distinct sources**, thus allowing for the formation of the most comprehensive expression data up to date.
在xena上面可以看到,TCGA和GTex、Target数据库的测序数据已经被重新计算整合到了一起,提供了各种格式的文件。
这里上游分析使用的是RSEM,而不是featurecout,导致得到的数据并不是标准的count值,是非整数。
1.RSEM 对应的差异分析包是EBSeq
RSEM (RNA-Seq by Expectation-Maximization)
关于它的下游分析,官网建议使用的R包是EBSeq:
EBSeq:http://www.bioconductor.org/packages/devel/bioc/html/EBSeq.html
但市面上公认最好的差异分析R包是DESeq2,limma,edgeR。有没有办法将RSEM的counts矩阵交给三大R包来处理呢?
2.能不能用limma和edgeR
这个问题刚好是关于TCGA的RSEM数据,
https://support.bioconductor.org/p/63981/#64004
limma的作者亲自回答了:
The RSEM expected counts from the TCGA project will work fine with either limma-voom or edgeR. However, with such a large number of samples, limma-voom is easily the best choice from a computational point of view.
limma-voom是更好的选择。
关于expected_count和norm_count在这里也有讨论,即edgeR只能用expected,vomm理论上可以使用norm_count(只是可以不是必须)。
3.能不能用Deseq2
https://support.bioconductor.org/p/94003/#94028
作者说最好的办法是用tximport进行转换,其次就是四舍五入对矩阵进行取整,然后用 DESeqDataSetFromMatrix去导入即可。
4.tximport是Deseq2作者写的R包
tximport:将各种上游分析软件得到的数据转换给三大R包使用。
其中讲了如何将ERR格式的数据导入R,并生成矩阵。DESeq2 和edgeR除了需要count矩阵,还需要一个length矩阵,而limma则是需要经过一步 "scaledTPM" 或"lengthScaledTPM"转换,需要另外一个矩阵来做校正。在xena中的数据我们只能拿到一个count矩阵,因此这个暂时用不上,但不妨碍它是个好东西。
5.胰腺癌的RSEM数据TCGA+GTEx联合分析
用三大R包差异分析分别做出的火山图和热图:
他们的结果取交集:
这个代码我已经上传到了github,https://github.com/xjsun1221/RSEM_with_limma_edgeR_Deseq2
。因为是根据作者的回复自己摸索的,受目前水平影响不能保证完全正确,以后也可能会进行修改,不建议作为标准答案来学习,供参考,如果发现有问题请发邮件到xjsun1221@163.com来反馈,提前感谢啦!我将在简书中更新。
另外github上有另外两人贡献的代码:
RSEM对接limma:https://github.com/NabilaRahman/RNA-Seq-Pipeline
RSEM对接DESeq2:https://github.com/yh154/rnaseq-rsem-star-deseq2-gsea