基因富集分析 是在一组基因中找到具有一定基因功能特征和生物过程的基因集，在研究差异表达基因、筛选基因的后续分析中经常使用。
基因集，也叫gene set，也就是一系列具有相同功能的基因构成的集合，比如某一条代谢通路（pathway），其中有很多的基因，因此位于同一条通路下的基因就构成了一个基因集合。
A: 基因功能富集分析输入的数据来源包括基因列表及基因表达谱;
B: 基因功能富集分析输入的基因功能注释数据库包括GO, KEGG, MSigDB等;
C: 过代表分析ORA, 先计数基因列表与基因功能集共同的基因, 并利用2×2 列联表, 进行Fisher 精确检验;
D: 功能集打分FCS,以GSEA 为例, 从基因表达谱中对所有基因按照与表型的差异表达程度排序, 再利用KS 检验计算对待测基因功能集打分, 利用随机重排得到背景分布, 找出显著富集的基因功能集;
E: 基于通路拓扑结构PT, 以Pathway-Express 为例, 根据上下游关系对通路进行打分, 上游通路影响因子较大, 对待测基因功能集根据影响因子打分, 然后根据统计检验模型计算显著程度;
F: 基于网络的拓扑结构NT, 以NEA 为例, 将基因列表或基因表达谱和待测基因功能集放入生物网络中, 计算感兴趣的基因与待测基因功能集在网络之间的最短距离, 通过对网络重调得到背景分布, 找出显著富集的基因功能集;
G: 基因功能富集分析输出的结果是排序后显著富集的基因功能集, 并对相似度较高的基因功能集进行去冗余
基因本体数据库 (Gene Ontology)，旨在建立基因及其产物知识的标准词汇体系，涵盖了基因的细胞组分 (cellular component)、分子功能 (molecular function)、生物学过程(biological process)。Term是GO里面的基本描述单元。
- 分子功能(Molecular Function)：描述在个体分子生物学上的活性，如催化活性或结合活性。
- 生物学过程(Biological Process)：由分子功能有序地组成的，具有多个步骤的一个过程。
- 细胞组件(Cellular Component)：指基因产物位于何种细胞器或基因产物组中(如糙面内质网，核糖体，蛋白酶体等)，即基因产物在什么地方起作用。
GO注释系统是一个有向无环图 (Directed Acyclic Graphs)，包含三个分支(即细胞组分、分子功能 、生物学过程)，（有向指的是term之间的单向指向性关系，比如termA是内质网，termB是细胞器，规定A是B，却不能说B是A；无环指的是从任何一点开始沿着规定的指向都不能回到原点）。一般这种关系包括 is_a/part_of/regulates
Kyoto Encyclopedia of Genes and Genomes: 系统分析基因产物和化合物在细胞中的代谢途径以及这些基因产物的功能的数据库。
整合了基因组、化学分子和生化系统等方面的数据,包括代谢通路（KEGG PATHWAY）、药物（KEGG DRUG）、疾病（KEGG DISEASE）、功能模型（KEGG MODULE）、基因序列（KEGG GENES）及基因组（KEGG GENOME）等。KO是蛋白质或酶的一个分类体系，将同一条通路上功能相似、序列相似的蛋白质归为一类。
命名规则：
K（大写)+num 基因ID号，表示所有同源物种中具有相似结构或功能的一类同源蛋白，如：K04456表示丝氨酸蛋白激酶
ko+num 代谢通路，表示特定的生物路径，如：ko04151表示PI3K-Akt信号通路【也是我们常用的代谢通路】
M+num 表示模块，如：M00676表示PI3K-Akt信号模块
C+num 表示化合物，如：C00533表示一氧化碳
ECx.x.x.x 表示酶，如：EC2.7.11.1表示丝氨酸
R+num 表示反应名称
Over-representation Analysis 过表达分析,又称为"2X2法"
具体步骤：
-
获得一组感兴趣的基因(一般是差异表达基因)。
-
给定的基因列表与某个通路中的基因集做交集，找出其中共同的基因并进行计数(统计值)。
-
利用统计检验的方式来评估观察的计数值是否显著高于随机，即待测功能集在基因列表中是否显著富集。
(最常用的统计检验包括：超几何分布、卡方检验、二项分布。)
KEGG通路hsa04310指的是Wnt signaling pathway，包括160个基因，背景基因使用了8017个**；得到的差异基因数是825个，在KEGG能注释上的有337个，其中就有19个是这个通路的，概率高达5%（enrichKEGG方法都是用能在KEGG注释上的基因，比如这里用80而不是547），那么这个通路是不是在下调基因中被显著改变？超几何分布检验，得到p值。hsa04310通路在背景基因中查到的概率显著低于差异表达基因的。
优点：**
基于完备的统计学理论, 结果稳健、可靠。
缺点：
(1)仅使用了基因数目信息，而没有利用基因表达水平或表达差异值，为了获得感兴趣或者差异表达基因，需要人为的设置阈值；
(2)ORA法通常仅使用最显著的基因，而忽略差异不显著的基因。在获得感兴趣的基因时, 往往需要选取合适的阈值, 有可能会丢失显著性较低但比较关键的基因, 导致检测灵敏性的降低；
(3)将基因同等对待，ORA法假设每个基因都是独立的，忽视了基因在通路内部生物学意义的不同(如调控和被调控基因的不同)及基因间复杂的相互作用；
(4)ORA假设通路与通路间是独立的，但这个前提假设是错误的。
Functional Class Scoring 功能集打分
具体步骤：
1.根据案例和对照状态下的基因表达谱对基因组中所有基因表达水平的差异值进行打分或排序，或直接输入排序好的基因表达谱。
2.把待测基因功能集中的每个基因的分数通过特定的统计模型转换为待测基因功能集的分数或统计值。
3.利用随机抽样获得的待测基因功能集统计值的背景分布来检验实际观测的统计值的显著水平，并判断待测基因功能集在案例和对照实验状态下是否发生了统计上的显著变化。
优点：
相较于ORA 法在理论上有明显突破, 考虑到了基因表达值的属性信息, 以待测基因功能集为对象来进行检验, 也使得检验结果更加灵敏。
缺点：
(1)仍独立分析每一条通路，但同一个基因可能涉及多条通路，所以不同通路间的基因出现重叠，所以别的通路可能由于重叠的基因，也出现显著富集；
(2)仍然把待测基因功能集中的每个基因作为独立的个体, 忽略了基因的生物学属性和基因间的复杂相互作用关系。
Pathway Topology 通路拓扑结构
把基因在通路中的位置(上下游关系)，与其他基因的连接度和调控作用类型等信息综合在一起来评估每个基因对通路的贡献并给予相应的权重，然后再把基因的权重整合入功能富集分析。不同的PT方法在具体的权重打分时，采用了不同的方式。
GO 等注释数据库中基因功能集中不包含任何拓扑结构信息，仅提供了可能属于同一通路的所有基因列表。所以，PT方法不能被用于GO通路的富集分析。
优点：
对于研究较完善、拓扑结构完整的通路，基于PT的基因功能富集算法会有更强的显著性。
缺点：
对于通路拓扑结构存在依赖性，对于研究较少、信息不完善的通路稳健性较差，目前通路注释的不完善。
Network Topology 网络拓扑结构
目前NT法有一些不同的思路：
1.基于生物网络拓扑结构的富集分析方法，利用数据库中的基因相互作用关系来间接地把基因的生物学属性整合入功能的富集分析。
2.利用网络拓扑结构来计算基因对特定生物通路的重要性并给予相应的权重， 然后再利用传统的ORA 或 FCS 方法来评估特定生物通路的富集程度，如 GANPA 和 LEGO 等；
3.直接把基因列表中的功能富集问题利用网络转化为基因对的功能富集问题，如 NOA 等。
优点：
与传统方法相比，基于网络的基因功能富集分析方法加入了系统层面的基因重要性程度及关联信息，使得预测结果更加准确可靠。
缺点：
算法过于复杂，计算速度较慢。
Gene Set Enrichment Analysis (基因集富集分析)用来评估一个预先定义的基因集的基因在与表型相关度排序的基因表中的分布趋势，从而判断其对表型的贡献。其输入数据包含两部分，一是已知功能的基因集 (可以是GO注释、MsigDB的注释或其它符合格式的基因集定义)，一是表达矩阵 (也可以是排序好的列表)，软件会对基因根据其与表型的关联度(可以理解为表达值的变化)从大到小排序，然后判断基因集内每条注释下的基因是否富集于表型相关度排序后基因表的上部或下部，从而判断此基因集内基因的协同变化对表型变化的影响。
GO富集分析是先筛选差异基因，再判断差异基因在哪些注释的通路存在富集；这涉及到阈值的设定，存在一定主观性并且只能用于表达变化较大的基因，即我们定义的显著差异基因。而GSEA则不局限于差异基因，从基因集的富集角度出发，理论上更容易囊括细微但协调性的变化对生物通路的影响，尤其是差异倍数不太大的基因集
给定一个排序的基因表
L和一个预先定义的基因集S(比如编码某个代谢通路的产物的基因, 基因组上物理位置相近的基因，或同一GO注释下的基因)，GSEA的目的是判断S里面的成员s在L里面是随机分布还是主要聚集在L的顶部或底部。这些基因排序的依据是其在不同表型状态下的表达差异，若研究的基因集S的成员显著聚集在L的顶部或底部，则说明此基因集成员对表型的差异有贡献，也是我们关注的基因集。
GSEA计算中几个关键概念：
-
计算富集得分 (ES, enrichment score). ES反应基因集成员s在排序列表L的两端富集的程度。计算方式是，从基因集L的第一个基因开始，计算一个累计统计值。当遇到一个落在s里面的基因，则增加统计值。遇到一个不在s里面的基因，则降低统计值。
每一步统计值增加或减少的幅度与基因的表达变化程度（更严格的是与基因和表型的关联度，可能是fold-change，也可能是pearson corelation值）是相关的，可以是线性相关，也可以是指数相关。富集得分ES最后定义为最大的峰值。正值ES表示基因集在列表的顶部富集，负值ES表示基因集在列表的底部富集。
-
评估富集得分(ES)的显著性。通过基于表型而不改变基因之间关系的排列检验 (permutation test)计算观察到的富集得分(ES)出现的可能性。若样品量少，也可基于基因集做排列检验 (permutation test)，计算p-value。
-
多重假设检验校正。首先对每个基因子集s计算得到的ES根据基因集的大小进行标准化得到Normalized Enrichment Score (NES)。随后针对NES计算假阳性率。（计算NES也有另外一种方法，是计算出的ES除以排列检验得到的所有ES的平均值）
-
Leading-edge subset，对富集得分贡献最大的基因成员。
结果解读：
-
ES是Enrichment score 富集得分：表示基因集的基因在基因列表两端的富集程度。开始从左到右基因集的第一个基因开始，每次计算统计值并且逐次累加，有一个在基因列表中的基因，就增加计算的统计值，反之就减小。并且每次的增加和减小都被记录在上图红线的位置，它代表了基因和表型的相关性。结果得到的峰值计作富集得分ES，正值表示基因集的基因主要在基因列表顶部富集，负值表示在底部富集
-
绿色圆圈表示的Leading edge subset：表示对ES富集得分最大的基因集中的基因，ES为正值表示是峰左侧的基因，负值表示右侧的基因。这个图中的每一条线表示基因集中一个基因
第一步 基因排序：**
GSEA分析的第一步就是利用所有基因的表达数据，然后计算每个基因在两个分组（或者表型）ClassA、ClassB中的差异度（GSEA提供了6种算法，默认方法是signal2 noise，GSEA官网有提供公式），然后按照在两个表型的差异度从大到小排序，形成一个排好序的基因列表。
第二步 分析基因集是否富集：
这里的基因集，是事先根据功能或者其他的一些原理把很多的基因分类成不同的基因集合，具体一个基因集可以是某一个通路或者go term中的所有基因，也可以是一个miRNA靶标对应的多个基因。GSEA提供了多个分类基因集，在分析数据时，只需要选择不同基因集就可以，当然也可以自己制作基因集。我们可以对每一个小的基因集（GeneSet ）里面的基因对应一下上一步排序表里面的位置。
第三步：
就是计算基因集的ES值（enrichment score），之后对基因集的ES值进行显著性检验及多重假设检验，从而计算出显著富集的基因集。
与常规基于超几何检验的基因功能富集分析方法相比，GSEA分析有如下的优点：
-
不需要对基因进行差异显著的筛选，这样能保留那些表达变化不大，但是功能重要的基因，而传统的GO和KEGG富集分析是针对有差异的基因进行富集分析，相比之下GSEA分析保留了更多信息。举个例子：这里我们找到了两个基因在细胞增值通路里有显著差异，同时呢，也找到了两个基因在细胞凋亡通路里有显著差异，这时候用传统的方法无法确认我们研究的细胞表型是和凋亡还是和增值相关？这时候GSEA分析的优势来了，我们用GESA进行富集分析发现凋亡通路里面的基因除了那两个显著差异的基因表达显著上升外其他的相关基因也有整体的上升，而增值通路里面的基因就没有这种现象，因此我们就可以确定我们研究的表型是和凋亡相关的。
-
分析的是基因集而不是单个的基因，因为生物体要出现表型差异，要找到与表型差异相关的基因，单单通过差异分析是不够的，有时候甚至得到假阳性的结果，因为生物体出现某种表型（一两个基因表达存在差异）往往会有一系列与之相关的上游或者下游的基因发生变化，但不一定会有显著差异，因此我们对功能相关的基因作为一个整体做GSEA分析，比较集合中基因整体的表达量差异得到的分析结果更可靠。
-
目前GSEA提供的功能基因分类数据库有以下8种，主要是与人类基因相关的分类数据库MSigDB，动植物目前没有，所以GSEA的分析方法大多在人类相关的研究中应用。但是如果你可以自行按照GSEA官方说明制作基因功能分类数据库，就可以应用到任何动植物了。
seqtk sample harm_pep.fa 100 > test.fa