写在前面
时间拨回去 2015年,那时我接触生信已有一年,TBtools 开发尚在萌芽阶段。那会,我写了几款小的软件,包括 “blast3go”,为的是应对即将收费的“blast2go”。当然,后来相关功能都整合到 TBtools 中。而其中有一个重点功能,即 GO 富集分析。那会在 Bioinformatics中国 社群,我们开始了理论上是国内最早的公开社群学术Seminar(网络直播),我在其上也分享了相关学习经验和知识。主要介绍了最常用的基于超几何分布的富集分析方法,事实上,可以理解为简单的 Fisher's Extact Test。同时,也介绍了 GSEA,全称 Gene Set Enrichment Analysis。这个有趣。因为对于新手来说,概念似乎有点混乱。
常见的富集分析逻辑
目前,在科研工作中最为常见的三种富集分析策略,分别是:
- ORA,全称,Over Represent Analysis,即前述基于超几何分布的统计检验。TBtools当前使用的也是,这个方法,最大的优点是,非常快!可以说,适用于几乎所有情况,毕竟生物学问题,我们常常希望的是抓主要矛盾,而我们能做的,也是抓主要矛盾。
- 卡方检验,这个比较简单,针对任何一个GO条目或者KEGG通路,更或者其他注释,很容易拿到背景属于和不属于该条目的两个数字,选择集属于和不属于该条目的两个数字,于是就是列联表检验。这个方法其实用得不过,尽管也可以看到。
- GSEA,翻译,基因集合富集分析。Emmm,这个有点意思。一般来说,我们在现实生活中总是说 GO富集分析,KEGG富集分析,统称基因集合富集分析。而GSEA这个方法本身翻译过来,也是基因集合富集分析。但其中指代区别,还是要知悉。
既生ORA,何生GSEA?
事实上,关于富集分析方法的选择和优化,是一个持久不衰的热门话题,包括近期更新的国人佳作 KOBAS-i 也提出了新的更完善的方法,感兴趣的朋友,建议了解。那为啥如此?原因不用说,肯定是旧方法还是有一些值得改进的地方嘛。其中ORA最明显的问题,就是“硬阈值”。ORA要求的输入是,背景注释+选择集。而选择集要么就是差异表达基因,要么可能就是GWAS关联分析出来的基因。具体到差异表达基因上,往往要求 pvalue必须低于 0.05?或者0.01?同时差异表达倍数要达到 2倍以上。这是比较有趣的,因为我们都知道 1.99999.... 倍的差异其实也算差异。这是 GSEA 方法解决的主要问题之一。
对于 GSEA 方法,一般我们可以提供 背景注释+基因排序信息表。后者或者物种所有基因的表达矩阵,由软件计算差异表达倍数或其他来排序以及作为权重;或者是直接输入物种所有基因在某个两两比较的差异表达倍数,Log2FoldChange值;其实还可以是其他信号强度数值。于是不存在“硬阈值”。解决了ORA的问题。同时其实还解决了另外一个问题,这里埋伏笔,下述给出。
GSEA~真香定律
事实上,15年我用过两三次,做完一个医学项目之后,就再也没碰过 GSEA,因为我几乎在所有项目中发现,ORA的表现已经足够好,根本不需要GSEA。尤其是 GSEA 计算起来,真的非常非常非常慢(当然,非常感兴趣的可以了解一下大湿兄的R语言包 ClusterProfiler,据说很快)。由于我确实不懂 R语言,所以一直也没去折腾,索性就一直没用。
但最近做到一组数据,做的是果皮软化。很明显,果皮肯定是软化的,表型非常明显,细胞壁代谢相关基因必然也差异表达,但是在 常规富集分析中,没有任何相关显著富集结果。
这个不应该出现。我回头看了下差异表达基因列表,发现其实倍数差异排前的,几乎都与细胞壁代谢相关。直觉告诉我,果实软化过程,基因表达变化过大(对应差异表达基因数目偏多),导致没有合理的富集结果。于是,我想起来 GSEA。前述提及,他解决“硬阈值”的问题,但这个问题在本例其实没啥用,因为都差异表达的。而他其实同时还解决了一个问题,那么就是权重!只要我们使用每个基因的差异表达倍数作为权重,那么差异倍数大的,自然权重就高,可以被正确反映到富集分析结果上。
打包 GSEA 成 TBtools 插件
既然 GSEA 还是能用到我的日常项目中,我简单想了下,那不如~~~ 还是真香定律一把。索性打包成插件,相信也会有其他朋友跟我一样遇到类似的问题。既然打包成插件,那么使用就应该足够简单。TBtools 目前使用的 GO 富集分析功能,界面和使用方法如下:
而对于 GSEA 插件,界面和使用方法几乎一样,只需要,调整一个输入文件:
选择集
->所有基因排序信息
。大体如下,两列,第一列是基因ID ,第二列是排序信息,可以直接整理出来所有基因的 Log2FoldChange值,也可以是其他信号值。注意,必须带表头,如下
界面和使用方法如下:
同样,点击 Start 就可以了。
GSEA 结果简单解读
首先,GSEA 真的有点慢,建议跑上之后,玩 TBtools 的俄罗斯方块或者贪吃蛇,大概十几分钟之后就会有弹窗告诉你跑完了。
打开输出目录,可以看到非常多的文件。
打开网页可以看到
在 SnapSpot 中可以看到还不错的图片
但更多时候,我们会直接进入详情信息表
对于感兴趣的条目,可以点击进去查看详情。
GSEA使用的是置换检验,Permutation,感兴趣的建议了解下 PNAS 文稿。从这个结果来看,对该条目贡献最大的,即
Core Enrichment
中标记为 Yes
的基因。理论上来说,这些基因,即下游实验最值得关注的基因。
写在最后
我已经打包成 TBtools 插件,也自认为绝大多数人应该可以毫无难度是使用得了这个插件。但这个插件的只是个把小时的简单工作,暂时没有优化插件的计划,当然也要看到底大伙感兴趣的人有多少。于是,有两个点:
- 插件的安装,请直接在 TBtools 中 Install Plugin 菜单选择插件并安装;
- 插件的获取:暂时不放到插件商店。。。。