下面是2021六月份学徒Sophie的投稿
最近参加了生信技能树曾老师的GEO数据挖掘月学徒培养,对一些文章中的GSE数据集走标准化分析流程。
小洁老师在去除批次效应的探索文件里给出了两种方法,一个是用R包limma中的函数removeBatchEffect(),另一个是R包SVA中的函数ComBat()。
阅读相关文献查到,也有很多文章用了RobustRankAggreg找个R包进行不同芯片数据的差异基因整合。区别就在于,是先进行样本整合,然后去除批次效应,最后进行差异分析;还是,对每个数据集进行独立分析,然后RRA整合DEG,拿到多个数据集共有的DEG,所以我尝试比较了接触到的这两种方法的结果差异。
首先
我们应该先了解,这个多个芯片数据差异基因整合R包:RobustRankAggreg,介绍在:https://cran.r-project.org/web/packages/RobustRankAggreg/index.html。
然后
limma中的函数removeBatchEffect()和SVA中的函数ComBat()介绍很多,直接浏览器搜索关键词即可。
接下来就是结果啦~
- 对每个数据集进行独立分析后,用RRA整合,然后拿top10 上下调基因与文章中给出top10进行对比,结果如下(带星标的为一致的基因):
- 先整合所有数据集,然后用limma去除批次效应,再进行差异分析,然后拿top10 上下调基因与文章中给出top10进行对比,结果如下(带星标的为一致的基因):
这么看起来,似乎是 我们的用limma去除批次效应,再进行差异分析,然后拿top10 上下调基因与文章的结果一致性好一点哦。
其实是因为从数量上来说,用RRA整合得到的差异基因数量为39个(数据集命名为list_RRA),用limma整合得到的差异基因数量为858个(数据集命名为list_limma),如果一个方法拿到的基因数量足够多,当然是更有可能保护文献的基因集。
那两个差异基因的关系是怎样的呢?见下图
所以结论是,RobustRankAggreg处理,拿到的DEG少,因为没有进行样本的整合,获得DEG所用的样本数量少。(换句话说,是比较严格的挑选了上下调基因,是在多个数据集都很保守的上下调差异基因)
但是优点在于,二分组样本的差异分析流程,对数据前期处理是最简单的,只需要根据pd对样本进行合理分组,即可。而对于limma包的函数removeBatchEffect()来说,虽然前期处理会麻烦一点,但是,能拿到更多的差异基因,因为扩大了样本量。而且,就目前这篇文章的数据集分析结果来说,与文章数据分析结果一致性更高。
还有一个问题
那到底是拿到多的DEG好,还是少的DEG好呢?其实关键不在多和少,而在于根据我们的背景知识,能不能找到符合预期的靶基因,如果RRA不行,那就换一种方法呗~