详解:基因集富集分析GSEA

1) Why GSEA?

传统KEGG(通路富集分析)和GO(功能富集)分析时,如果富集到的同一通路下,既有上调差异基因,也有下调差异基因,那么这条通路总体的表现形式究竟是怎样?是被抑制还是激活?或者更直观点说,这条通路下的基因表达水平在实验处理后是上升了呢,还是下降了呢?

传统的富集分析,针对总体的差异基因,不区分哪些差异基因是上调还是下调。

这是因为,传统的富集分析根本不考虑基因表达量的变化趋势,其算法的核心只关注这些差异基因的分布是否和随机抽样得到的分布一致,即使后期绘图时,在通路图上用不同颜色标记了上下调的基因,但,由于没有采用有效的统计学手段去分析某条通路下的差异基因的总体变化趋势,这使得传统的富集分析结果无法回答上述的问题。

于是,有人提出,分别提取上调或者下调的差异基因来进行传统富集分析,由于事先根据上调下调(表达量变化趋势)对差异基因进行了筛选,从而回避了上面的问题。但,这样的做法有失偏颇,因为费舍尔精确检验就是想要证明这个差异基因列表不是随机抽样得到的,而事先对差异基因列表的进行上调或下调的过滤,就会对结果的随机性造成了干扰,最后得出的结论其准确性也受影响。而且,在细胞内发挥生物学功能时,上调和下调基因是共同发挥作用,进行富集分析时,将上下调分开进行分析,也不符合实际情况。

想象一下,上调基因和下调基因分开富集,然后富集到了同一条通路,这怎么解释?

所以,准确来说,传统的富集分析只能定位到功能,这些差异基因与哪些功能相关,而不能回答某条通路的整体表现形式。想要回答这个问题,所以我们需要GSEA富集方法的结果

2)What is GSEA?

GSEA(Gene Set Enrichment Analysis):基因集富集分析,由Broad Institute研究所提出的一种富集方法,同时还提供对应的分析软件GSEA和一个基因集数据库MSigdb[http://software.broadinstitute.org/gsea/msigdb/index.jsp]

GSEA的原理图

GSEA的输入是一个基因表达量矩阵,其中的样本分成了A和B两组,首先对所有基因进行排序,简单理解就是根据处理后的差异倍数值进行从大到小排序, 用来表示基因在两组间的表达量变化趋势。排序之后的基因列表其顶部可看做是上调的差异基因,其底部是下调的差异基因。

GSEA分析的是一个基因集下的所有基因是富集在这个排序列表的顶部还是底部,如果在顶部富集,可以说,从总体上看,该基因集是上调趋势,反之,如果在底部富集,则是下调趋势。

GSEA分析有三个特点:

1. 分析的基因集合而不是单个基因;

2. 将基因与预定义的基因集进行比较;

3. 富集分析;

假设1个比较组“MUT vs WT”的差异gene集(MUT为实验组,WT为对照组),进行GSEA富集分析,结果如下图:

GS:基因集(通路)的名字。

SIZE:代表该基因集(通路)下的基因总数。

ES:代表Enrichment score,NES代表归一化后的Enrichment score。

NOM p-val:代表p值,表征富集结果的可信度。

FDR q-val`代表q值, 是多重假设检验矫正后的p值,注意GSEA采用pvalue < 5%, qvalue < 25% 对结果进行过滤。

对于某个基因集下(通路里)的每个基因给出了详细的统计信息,如下图

RANK IN GENE LIST:代表该基因在排序中的位置。

RANK METRIC SCORE:代表该基因排序量的值,即:处理后的foldchange值。

RUNNIG ES:代表累计的Enrichment score。

CORE ENRICHMENT:代表是否属于核心基因,即对该基因集的Enerchment score做出了主要贡献的基因。

上图表格中的数据对应下面这张图

Enrichment Score的折线图

上图分为3部分,如下:

第一部分:最顶部的绿色折线为基因Enrichment Score的折线图。纵轴为对应的Running ES, 在折线图中有个峰值,该峰值就是这个基因集的Enrichemnt score,峰值之前的基因就是该基因集下的核心基因。横轴代表此基因集下的每个基因,对应第二部分类似条形码的竖线。

第二部分:类似条形码的部分,为Hits,每条竖线对应该基因集下的一个基因。

第三部分:为所有基因的rank值分布图,纵坐标为ranked list metric,即该基因排序量的值,可理解为“公式化处理后的foldchange值”。

我们假设是针对比较组“MUT vs WT”的差异gene集进行分析,MUT为实验组,WT为对照组,差异gene的差异倍数计算公式为:

通常统计时,对foldchange取log值(取对数)。如果log2(foldchange)>0,表明实验组表达量高于对照组,即,实验组相对于对照组上调。如果log2(foldchange)<0,表明实验组相对于对照组下调。
上图“Enrichment Score的折线图”对应的纵轴值全部大于0,显示这个基因集是在MUT组高表达的(即,此基因集(通路)在MUT组上调),下图是一个在WT组中高表达的示例。

下图中其Enrichment score值全部为负数,其峰值右侧的基因为该基因集下的核心基因。

对于一个基因集而言,定义其中对Enrichment score贡献最大的基因为核心基因,也称之为leading edge subset, 参考下图

对于Enrichment score为正数的基因集而言,其核心基因是峰值之前的基因,对于Enrichment score为负数的基因集而言,其核心基因是峰值之后的基因。

上表KEGG富集分析结果,最后一列为LEADING EDGE, 在这一栏中,包含以下3个统计量:tags、list、signal
tags表示核心基因在该基因集基因总数的占比,而list表示核心基因占所有基因总数的比例,signal利用这两个指标计算得到,公式如下

N代表所有基因的数目,Nh代表该基因集下的基因总数。对于一个基因集而言,当核心基因的数目和该基因集下的基因总数相同,signal取值最大,当该基因集的基因数目和所有基因数目接近时,signal的取值接近于0。

GSEA软件默认的输入是基因集合(基因表达量矩阵和样本分组),然后内置的进行归一化,进行差异分析,计算singal2noise等统计量,其本质就是进行差异分析,计算出类似foldchange的统计值,但,GSEADE 归一化算法是否适用于我们输入的表达量矩阵,在计算基因的foldchange值时有没有考虑生物学重复本身的变化程度,这些都导致其计算出的foldchange值并不能满足我们的需求,更加有效的做法是采用专用的差异分析软件(比如DESeq2、Edger等)计算出的foldchange值来进行富集分析。

GSEA的核心是Enrichment score的计算,除了GSEA软件外,还有很多的工具也都支持这个算法,如果想要利用DESeq2等工具自定义计算出的基因排序列表进行富集分析,更推荐使用cluster Profiler等第三方工具。

GSEA富集过程包括三步骤:

1. 计算富集分数(Enrichment Score);

2. 估计富集分数的显著性水平;

3. 矫正多重假设检验;

GSEA的具体原理请见PNAS文章:Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A. 2005 Oct 25;102(43):15545-50. Epub 2005 Sep 30.

3)小结

常规富集分析侧重于比较两组间的基因表达差异,集中关注少数几个显著上调或下调的基因,这种方式可能由于筛选参数的不合理导致漏掉一些关键信息,比如,部分差异表达不显著却有重要生物学意义的基因。因为实际通过测序检测到的RNA表达变化,是层层的负反馈调控后的结果,并且不同组织对于表达差异的敏感度是不同的:在神经递质系统内,一个1.2倍的表达差异即可能产生极其显著的效应。

下图示例中,无论设置怎样的差异筛选条件,如果仅做普通的Pathway富集分析的话,一定会漏检至关重要的Myc通路。这个示例非常典型,不仅在于Myc作为重要的癌基因广为人知,并且这里Myc在实验条件下活性改变后引起的下游基因表达变化,非常具有代表性:即并非所有的下游基因都会展现出强烈的表达改变,但它们会呈现出一致的趋势。GSEA的优势在于,无需做差异分析,直接用所有基因的表达量进行分析,可检测出不显著但是一致的差异表达趋势的基因集。

GSEA也可用于判断,某条通路在某组样本中是激活还是抑制。

当然,常规富集分析和GSEA分析没有说哪个更好,实际应用中能解决问题即可,引用一句名言:黑猫白猫,能抓住老鼠的就是好猫。

--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------I am a line ! Thanks!----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------

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

推荐阅读更多精彩内容