学习基因通路富集分析软件GSEA

GSEA是一种无阈值方法,可根据其差异表达等级或其他分数对所有基因进行分析,无需事先进行基因过滤。当基因组中的所有或大多数基因(例如,RNA-seq数据)可获得rank时推荐使用GSEA进行通路富集分析, 然而,当仅有一小部分基因具有rank可用时,如,在确定显着突变的癌症基因的实验中,GSEA并不合适。GESA使用基于一个置换矩阵检验来分析rank gene list。GSEA搜索基因在rank gene list的顶部或底部富集通路,这比单凭偶然的机会所能预料到的还要多。 例如,如果最顶端的差异表达基因参与细胞周期,这表明细胞周期通路在实验中受到调控。相反,如果细胞周期基因在整个rank gene list中随机分散,则细胞周期途径可能不会受到显着调节。要计算通路的富集分数(ES),GSEA逐步从顶部到底部检测rank list的基因,如果基因是该通路的一部分则增加ES,否则降低分数。这些运行总和值是加权的,因此放大了top(和bottom)排序基因的富集,而中等水平基因的富集则没有被放大。ES分数被计算为运行总和的最大值并相对于通路大小进行归一化,从而得到标准化的富集分数(NES),其反映了列表中通路的富集。正NES值和负NES值分别表示列表顶部和底部的富集。最后一个基于置换的p值被计算,并用多次测试进行矫正以产生基于置换的错误发现率Q值,Q值的范围从0(非常显著)到1(不显著)。从排序基因列表的底部开始进行相同的分析,以鉴定在列表底部富集的通路。使用FDR Q值阈值(例如,Q <0.05)选择所得到的通路并使用NES进行rank。 此外,GSEA分析的“leading edge”方面确定了对检测到的通路富集信号最有贡献的特定基因。

GSEA有两种确定ES的统计学显着性(P值)的方法:基因集置换和表型置换。 基因集置换测试需要rank list,并且GSEA将观察到的通路的ES与通过用随机取样的匹配大小的基因集(例如,1,000次)重复分析而获得的分数分布进行比较。表型置换测试需要所有样品的表达数据(例如,生物学重复),以及被称为“表型”的样品组,该方法是彼此之间比较(例如,病例与对照;肿瘤与正常样品)。对于具有有限突变和生物学重复的研究(即每种条件2至5次),推荐使用基因集置换。在这种情况下,差异基因表达值应在GSEA之外计算,使用包括方差稳定性的方法(例如edgeR ,DESeq 和limma / voom ),并在通路分析之前导入GSEA软件,表型置换应该与使用更多次重复(例如,每种条件至少10次)。表型置换方法的主要优点在于与基因集置换方法相比,它在排列过程中保持了具有重要生物学意义的相关性基因的基因组结构。

input data :数据是被TCGA鉴定的两种卵巢癌两种亚型差异表达的gene list。该rank先前基于基因表达数据分层为四种分子亚型,定义为分化,免疫反应,间充质和增殖。GSEA需要具有基因分数的RNK文件,该rnk文件有两列,第一列是基因ID,第二列是基因分数。基因组中的所有(或大多数)基因需要具有分数,并且基因ID需要与GMT文件中使用的基因ID匹配。

Load Data:将需要分析的数据加载进来,同时也将进行通路分析的基因集(GMT)加载进来。


Load Data

在下图的这个地方就可以看到你加载进来的数据

已经加载进来的数据

点击左边工具栏的Run GSEAPreranked


GSEA自己也提供了基因集文件,可以直接通过MSigDB资源从GSEA端口直接访问,不需要输入GSEA中。要定义GMT文件,可以在Select one or more genesets dialog对话框的第一个选项卡Gene Matrix(from website)中找到MSigDB基因集文件。如下图,而我自己提供了GMT文件,就选择Gene matrix (local gmx/gmt)选项卡,这下面有你在Load data的时候加载进来的GMT文件。

Number of permutations: 这指定了基因集随机化以创建空分布以计算P值和FDR Q值的次数。 使用默认值1,000个排列。

更多的置换次数需要更长的计算时间。 为了计算每个gene set的FDR Q值,通过置换每个基因组中的基因并重新计算随机组的P值来随机化数据集,此参数指定完成此随机化的次数。执行的随机化越多,FDR Q值估计就越精确(达到极限,因为最终FDR Q值将稳定在实际值)。 

rank list: 单击最右侧的箭头并突出显示排名文件,选择排序的基因列表。

单击Basic fields 的Show可以展示出其他选项:


Basic fields

Analysis name: 默认是‘my_analysis',也可以修改为自己想要的。

Max size: exclude larger sets:默认情况下,GSEA将上限设置为500,将此值设置为200以从分析中删除较大的集合。


Max size: exclude larger sets设置为200

Save results in this folder:默认的结果保存路径是:/home/cmy/gsea_home/output/aug25,也就是会默认保存在你的home目录下。可将其改为自定义的路径


我将其改为自定义路径

除了上面Basic fields 外,下面还有一个Advanced fields,点击右边的show就可以展示出来,在这里可以进行一些高级的设置。在高级设置里面我没有修改,都是GSEA默认的。

Advance fields

设置结束之后点击下面"run”箭头按钮:


最下面右边有一个run的按钮

在run的过程中你会看到下图这个最右边的数字会变化,而且左下角的GSEA report 中Status 会变成running. run的时间长短和速度的快慢取决于你的电脑配置。

running


GSEA reports Status

运行结束之后GSEA report 中Status会从running 状态变成success,点击Success就可以以网页的形式查看你的结果,根据上调还是下调会展示为两组。结果文件会保存在Save results in this folder设置的文件夹中。


run结束


web展示的结果,上调为第一组结果,下调为第二组展示结果


设置的保存结果的文件夹下会多出两个结果文件夹

GSEA分析的结果通路富集排列在top的genes为上调的基因, na_pos(na表示'not available',因为我没有将表型标签的cle格式的文件输入,因此默认为na,而pos表示positive,neg:表示下调(negtive))


上调与下调结果概括

对于上图的结果解释可见下图:


上调组基因结果图的解释

点击enrichment results in html,可以在网页查看高表达基因集的富集的结果,如下图:


enrichment results in html结果展示

GS:基因集的名字,SIZE:基因集下的基因总数,ES:Enrichment score, NES:归一化后的Enrichment score,  NOM p-val:p-value,表征富集结果的可信度,FDR q-val:q-value, 是多重假设检验矫正后的p值,注意GSEA采用pvalue < 5%, qvalue < 25% 对结果进行过滤。

点击GS DESC可以跳转到每个基因集详细结果页面,如下所示:

GS DESC下的表

首先是一个汇总的结果,Upregulated in class说明该基因集在na-pos这组中高表达,其他的信息和上一副图中的表描述的一样,除此之外,还有详细的表格,如下所示:


GS DESC下的表

上表对于该基因集下的每个基因给出了详细的统计信息,RANK IN GENE LIST代表该基因在排序号的列表中的位置, RANK METRIC SCORE代表该基因排序量的值,比如foldchange值,RUNNIG ES代表累计的Enrichment score, CORE ENRICHMENT代表是否属于核心基因,即对该基因集的Enerchment score做出了主要贡献的基因。

该表格的基因集对应下面这张图:

GS DESC下的图

该图分为3个部分:

第一部分为基因Enrichment Score的折线图,横轴为该基因集下的每个基因,纵轴为对应的Running ES, 在折线图中有个峰值,该峰值就是这个基因集的Enrichemnt score,峰值之前的基因就是该基因集下的核心基因,即对该基因集的Enerchment score做出了主要贡献的基因。

第二部分为hit,用黑色线条标记位于该基因集下的基因

第三部分为所有基因的rank值分布图

从上图可以看到,其Enrichment score值全部为正数,对应的在其峰值左侧的基因为该基因集下的核心基因。

在总的html页面中,还给出了如下信息:

html页面结果下的一些描述信息

Dataset details给出了基因总数,Gene Set details给出了基因集的信息,默认根据基因集包含的基因个数是先对基因集进行过滤,最小15个,最大500个基因,(由于我自己设置了最大为200个),所以这里过滤掉了12032个基因集,剩余的 4714 个基因集用于分析。


学习与参考:

1.https://blog.csdn.net/weixin_43569478/article/details/83745105

2.NCBI - WWW Error Blocked Diagnostic

3.Reimand J, Isserlin R, Voisin V, et al. Pathway enrichment analysis and visualization of omics data using g:Profiler, GSEA, Cytoscape and

EnrichmentMap[J]. Nature Protocols, 2019, 14(2): 482-517.

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

推荐阅读更多精彩内容