关于GSEA富集结果条目数一跑一变的迷思和恍然大悟

当只有自己发生状况时,归结给运气和自己太菜;当别人发生同样的状况时,归结给共同的客观原因,比如网速太慢、电脑太烂;当伟大的互联网告诉你有更多的人遭遇了同样的困惑,才终于投向科学,怪罪无知,无知当然是自己的无知。

状况1: GSEA-GO 同样的设定返回了不同的结果

在做 jimmy 大神布置的作业:对比 Gorilla, clusterProfiler, topGO 三种工具

作业做到 GSEA 法

## clusterProfiler GSEA GO-BP
library(clusterProfiler)
library(hgu133plus2.db)
load("genelist_ENTREZID_Decr_pVal.Rdata")

### with decreasing p.val
gseaGO_BP <- gseGO(geneList     = entrezlist_dcrpv,
                   OrgDb        = hgu133plus2.db,
                   keyType      = "ENTREZID",
                   ont          = "BP",
                   nPerm        = 1000,   ## 排列数
                   minGSSize    = 5,
                   maxGSSize    = 500,
                   pvalueCutoff = 0.95,
                   verbose      = TRUE)  ## 不输出结果
gseaGO_BPresult <- gseaGO_BP@result
save(gseaGO_BPresult, file = "cp_gseaGO_BPresult.Rdata")

得到的结果是下面这样的,出现了65个条目和10个条目两种,似乎还有0条,但手慢无截图

当时的想法是这样的:

可能和网络也有关系?几次输入同样的参数,得到的结果并不一样....有时是报错 no term, 迷惑desu

毕竟日常网络不好,总天真的以为5G来了一切就解决了┑( ̄Д  ̄)┍

状况2: GSEA-KEGG 同样的设定返回了不同的结果

没有错,前后呼应,这是自称资质平平送老迎新的技能树资深半席讲师遇到的状况。依然是同样的参数,不同的结果。

讨论:也许这有一定的随机性?

后来果然发现 gseKEGG() 的一个逻辑参数 seed 是可以设置的,默认值为 FALSE, 那么就惊喜地设置为 TRUE, 又惊喜地发现了依然不同的结果。(sad

网友们的困惑

回到最初的问题,为什么每次会不一样呢?第二个问题,seed 这个参数到底有没有用?

自己思索无果就求助于互联网(然而其实应该别瞎想直接google

大部分第一个问题基本没有已解决的下文,第二个问题倒是终于发现了这位心明眼亮的选手:

当然也有人说并没有用....

不如自己动手尝试:

set.seed(1984)
gseaGO_BP3 <- gseGO(geneList     = entrezlist_FC,
                    OrgDb        = hgu133plus2.db,
                    keyType      = "ENTREZID",
                    ont          = "BP",
                    nPerm        = 1000,   ## 排列数
                    minGSSize    = 5,
                    maxGSSize    = 500,
                    pvalueCutoff = 0.05,
                    verbose      = TRUE,
                    seed         = TRUE)

第一次:no term enriched under specific pvalueCutoff...

第二次:结果又出现辽

看来seed确实是没用的。

另一位心明眼亮的网友如是说:

再回头去看,gseGO()gseKEGG() 默认的方法是 "fgsea".

We implemented GSEA algorithm proposed by Subramanian(Subramanian et al. 2005). Alexey Sergushichev implemented an algorithm for fast GSEA analysis in the fgsea(S., n.d.) package.

In DOSE(Yu et al. 2015), user can use GSEA algorithm implemented in DOSE or fgsea by specifying the parameter by="DOSE" or by="fgsea". By default, DOSEuse fgsea since it is much more fast.

https://yulab-smu.github.io/clusterProfiler-book/

而 fgsea 算法的文章里有几段话(反正公式是看不懂的也就看点人话猜一下这亚子

Gene set enrichment analysis is a very widely used method for analyzing gene expression data. It allows to select from an a priori defined list of gene sets those which have non-random behavior in a considered experiment.

The method has a major drawback of being relatively slow....That can be done in a straightforward manner by sampling random gene sets.

Instead of generating nm independent random gene set for each permutation and each gene set we will generate only n radom gene sets of size K.

The preranked gene set enrichment analysis takes as input two objects: an array of gene statistic values S and a list of query gene sets P. The goal of the analysis is to determine which of the gene sets from P has a non-random behavior.

Sergushichev A. An algorithm for fast preranked gene set enrichment analysis using cumulative statistic calculation[J]. BioRxiv, 2016: 060012.

大概就是说传统的GSEA方法很慢,可以通过对随机的基因集抽样解决这个问题,也就顾名思义是 a fast gene set enrichment analysis (FGSEA) 了。

与此同时四通八达的 jimmy 大神去问了神包 clusterProfiler 的作者,神通广大的Y叔——
Y叔回答:“nPerm次数多点就行 seed没用

毕竟👇

For each p ∈ P we need to find the enrichment statistic value and to calculate a p-value of this not to be random. To calculate a p-value for gene set p we can obtain an empirical null distribution by sampling n random gene sets of the same size as p.

Sergushichev A. An algorithm for fast preranked gene set enrichment analysis using cumulative statistic calculation[J]. BioRxiv, 2016: 060012.

所以,多试试总会得到想要的。


最后,向大家隆重推荐生信技能树的一系列干货!

  1. 生信技能树全球公益巡讲:https://mp.weixin.qq.com/s/E9ykuIbc-2Ja9HOY0bn_6g
  2. B站公益74小时生信工程师教学视频合辑:https://mp.weixin.qq.com/s/IyFK7l_WBAiUgqQi8O7Hxw
  3. 招学徒:https://mp.weixin.qq.com/s/KgbilzXnFjbKKunuw7NVfw
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 202,529评论 5 475
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 85,015评论 2 379
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 149,409评论 0 335
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,385评论 1 273
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,387评论 5 364
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,466评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,880评论 3 395
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,528评论 0 256
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,727评论 1 295
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,528评论 2 319
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,602评论 1 329
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,302评论 4 318
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,873评论 3 306
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,890评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,132评论 1 259
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 42,777评论 2 349
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,310评论 2 342