GSEA分析详细步骤

GSEA分组分析


这段时间有空为了某人学习了一下不属于我的领域的东西——GSEA分析(基因富文本分析),下面总结最近所学GSEA分析的步骤及方法。

1. 根据基因文件进行处理筛选

对于已有的基因文件,一般是xslx或者是csv格式(不要问基因文件哪来的,外行人的我猜测应该是可以从某某网站下载),而这个文件一般排序都是乱的且只有一列,所以这个时候就要进行分列并排序了。如图为基因文件原始格式:

在这里插入图片描述

接下来进行分列操作

  • 因为数据全部都在第一列,所以选中整个第一列
  • 点击数据-分列
  • 弹出分列窗口之后点击下一步,然后这里根据表格的形式选择对应的分隔符号,我这里选择Tab键就可以了,点击下一步,然后点击完成

接下来根据你需要分析的基因进行表达量大小的排序:

  • 选中你需要分析的基因那一列,这里我以HEATR1基因为例
  • 然后点击数据-排序
  • 确定选择扩展选定区域选项
  • 点击排序即可

接下来进行数据转置操作:

因为此时的数据行和列是反的,不利于我们进行GSEA分析,因此需要进行转置。

  • 新建一个空白表格
  • 复制刚刚分列的表格内的全部内容
  • 在空白表格以转置的格式粘贴

如图即为分列-排序-转置后的表格,瞬间清爽许多。

在这里插入图片描述

至此已经完成一半了,继续往下看。

2. 文件内容补充及改名

  • 因为步骤1所得到的文件GSEA软件并不能识别,为了能让GSEA软件识别,需要在表格第一行增添2空白行

  • 同样还需要在第二列增添一列空白列,完成行和列增添后的表格如图所示

    在这里插入图片描述

  • 接着按照如下格式在空白行、空白列中增添内容。这里解释一下12197代表的是表格的每行的基因个数(行基因数),30表示每列的样本个数(列样本数)。#1.2是固定格式。然后在第二列也就是30下面填入Description,Description下面填入NA,这也是固定格式。

    在这里插入图片描述

  • 接下来对以上文件进行保存,选择另存为,保存类型一定要选择文本文件(制表符分隔)(*.txt),文件名则切换输入法为英文,然后在文件名里面先输入双引号,双引号里面填你需要的名字并加上.gct。比如我甜的格式就是:

    "cell_heatr1_high_low_gct.gct"
    

    (命名的时候尽量都用下划线_而尽量少用中划线-,因为有些编程语言不识别中划线,而下划线基本都可以被识别。还有的人也许喜欢选择保存类型之后就随便的命名然后直接更改文件的后缀,这种方法大概率在后期文件导入到GSEA软件时会出现导入失败,尽量别用这种方法)
    在这里插入图片描述

  • 至此,GSEA所需的gct文件就已经完成了。接下来制作cls表型文件,制作格式如下图所示:

    在这里插入图片描述

    这里解释一下,在你对HEATR1基因的30个样本排序之后,你需要自己设置一个表达量的阈值,高于这个阈值的则为高表达,低于这个阈值的则为低表达,这也就是cls表型文件最下面的名称的由来

  • 最后另存为,保存类型一定要选择文本文件(制表符分隔)(*.txt),文件名则切换输入法为英文,然后在文件名里面先输入双引号,双引号里面填你需要的名字并加上.cls。比如我甜的格式就是:

    "cell_heatr1_high_low_cls.cls"
    

3. GSEA软件使用

  • 首先导入文件,点击左上角的Load data

  • 导入方式有三种,一般按照喜好自己选择,这里我使用方式1,点击Browse for files,选择之前制作好的gct文件和cls文件,导入成功如下图所示:

    在这里插入图片描述

  • 接着点击左上角的Run GSEA,在弹出的窗口中,首先选择表达数据集Expression dataset,可以看到只有我们导入的gct文件,所以选择它即可

  • 然后选择基因背景文件Gene sets dataset,我这里以kegg通路演示,具体根据你的需要来选择。

    在这里插入图片描述

  • 然后是置换检验的步数Number of permutations,默认选择1000,一般不进行修改,需要改的话根据你的需求自己改

  • 接下来选择表型文件Phenotype labels,表型文件直接选择刚刚导入的即可,第3步选择分组对比,这里随便选一个即可,比如我选了Low分组,那么我后期得到的正向通路就会和Low基因比较正相关。

在这里插入图片描述

  • 接下来选择是否要将基因转换为gene symbols格式,这里选择不转换

    在这里插入图片描述

  • 接下来选择置换检验的类型,一般样本量比较少时选择gene_set,该选项会影响结果的p值,所以选的时候注意一下
    在这里插入图片描述

  • 接下来是芯片平台Chip platform选择,因为我们不转换gene symbols文件,因此这里可以不要选择,如果要转化则要选择。该参数主要用来选择gene symbol与探针、或者gene ID相互对应的注释文件,当上边Collapse dataset to symbols参数设置不是No_Collapse时,此参数必须选择。这里,我们输入为gene symbol,因此不选

  • 最后是基本设置,具体设置见下图

    在这里插入图片描述

  • 最后点击最下面绿色的箭头run即可,最终结果就会在刚刚设置的文件路径下。如图所示

    在这里插入图片描述

GSEA单基因分析

单基因分析步骤和分组基因步骤不同的地方就两点,其余的步骤参考分组分析。

  • 基因表型文件Phenotype labels设置不同,具体操作见图片

    在这里插入图片描述

    在这里插入图片描述

    在这里插入图片描述

    至此,这里就设置好了

  • 接下来是第二处不同,也就是基因排序参数的修改,一般选择pearson,代表是基于其他基因与该基因的相关系数大小排序。(如有需要可自行更改其它排序方法)

    在这里插入图片描述

  • 接下来是第二处不同,也就是基因排序参数的修改,一般选择pearson,代表是基于其他基因与该基因的相关系数大小排序。(如有需要可自行更改其它排序方法)

  • 别的地方和分组分析都一样,设置完之后运行即可


https://www.bilibili.com/video/BV1CM4y1G7A3?from=search&seid=18177858094065419746&spm_id_from=333.337.0.0

<think>我们被要求生成关于GSEA富集分析的实现代码。根据引用内容,特别是引用[3]提供了具体的R代码示例,我们将以此为基础构建一个完整的GSEA分析流程。注意:用户要求的是GSEA富集分析的代码,因此我们将提供R语言的实现,使用clusterProfiler包。步骤:1.准备输入数据:需要一个排序的基因列表(按logFC等指标排序)和基因集(例如KEGG、GO或自定义基因集)。2.执行GSEA分析:使用GSEA函数。3.可视化:使用gseaplot2绘制富集结果。根据引用[3]的代码,我们注意到需要准备两个文件:term2gene(基因集与基因的对应关系)和term2name(基因集ID与名称的对应关系)。此外,输入基因列表需要按logFC从大到小排序。我们将提供两种情况的代码:情况1:使用内置的KEGG或GO数据库(如引用[1]所示,使用gseKEGG或gseGO函数)情况2:使用自定义基因集(如引用[3]所示,使用GSEA函数并传入TERM2GENE和TERM2NAME)由于用户没有明确说明使用哪种基因集,我们将分别展示两种情况。注意:引用[4]提到,分析结果需要关注p.adjust、NES和Description。代码将包括:-数据准备-执行GSEA-结果可视化-结果保存我们将按照以下结构组织代码:1.情况1:使用KEGG数据库2.情况2:使用自定义基因集注意:在情况1中,我们不需要自己构建term2gene和term2name,因为gseKEGG函数内部会处理。另外,引用[3]中使用了GSEA函数,并设置了nPerm=1000,pvalueCutoff=1(表示不根据p值过滤,后续可以自己过滤),seed设置可重复。开始编写代码:</think>以下是使用R语言进行GSEA富集分析的完整代码实现,基于clusterProfiler包: ```r # 加载必需的包 library(clusterProfiler) library(enrichplot) # 1. 准备输入数据 ----------------------------------------------------------- # 读取差异表达分析结果文件(需包含基因名和logFC值) data <- read.csv("diff_genes.csv", header=TRUE, sep=",") # 检查列名并标准化(确保第一列是基因名,第二列是logFC) colnames(data)[1:2] <- c("gene", "logFC") # 创建排序的基因列表(按logFC从大到小排序) genelist <- data$logFC names(genelist) <- data$gene genelist <- sort(genelist, decreasing=TRUE) genelist <- genelist[genelist != 0] # 移除logFC=0的基因 # 2. 选择分析模式 ---------------------------------------------------------- ### 模式1:使用内置数据库(KEGG或GO)### # KEGG分析(人类示例) gsea_result <- gseKEGG( geneList = genelist, organism = "hsa", # 人类:hsa,小鼠:mmu pvalueCutoff = 0.25, # 实际过滤padj值 pAdjustMethod = "BH", # Benjamini-Hochberg校正 seed = TRUE # 确保结果可重复 ) ### 模式2:使用自定义基因集 ### # 读取自定义基因集文件(两列:通路ID, 基因名) term2gene <- read.csv("custom_geneset.csv", header=TRUE)[, 1:2] colnames(term2gene) <- c("term", "gene") # 读取通路描述文件(两列:通路ID, 通路描述) term2name <- read.csv("pathway_names.csv", header=TRUE)[, 1:2] colnames(term2name) <- c("term", "name") # 执行自定义GSEA分析 set.seed(123) # 固定随机种子 gsea_result <- GSEA( geneList = genelist, nPerm = 1000, # 排列次数 pvalueCutoff = 0.05, TERM2GENE = term2gene, TERM2NAME = term2name, seed = TRUE ) # 3. 结果可视化 ----------------------------------------------------------- # 查看显著富集通路(按NES排序) head(gsea_result@result[order(gsea_result$NES, decreasing=TRUE), ]) # 绘制Top10通路GSEAgseaplot2(gsea_result, geneSetID = 1:10, # 前10个通路 title = "Top Enriched Pathways", pvalue_table = TRUE) # 显示P值表格 # 保存结果到文件 write.csv(gsea_result@result, "gsea_results.csv", row.names=FALSE) ``` ### 关键参数说明: 1. **输入文件要求**: - 差异基因文件:CSV格式,至少包含`gene`(基因名)和`logFC`(对数倍数变化)两列 - 自定义基因集:`term2gene`(通路ID与基因对应表),`term2name`(通路ID与名称对应表) 2. **核心参数**: ```r pvalueCutoff = 0.05 # 显著性阈值(推荐0.05) nPerm = 1000 # 排列次数(通常1000-10000) pAdjustMethod = "BH" # 多重检验校正方法 ``` 3. **结果解读要点**[^4]: - `NES`(标准化富集分数):|NES|>1表示显著富集,正负值分别代表上/下调 - `p.adjust`:校正后的P值(<0.05显著) - `core_enrichment`:核心富集基因列表 ### 注意事项: 1. 基因名需统一格式(如SYMBOL或ENTREZID) 2. 对于物种数据库(KEGG/GO),正确设置`organism`参数: - 人类: `"hsa"`,小鼠: `"mmu"`,大鼠: `"rno"` 3. 推荐使用**云分析平台**简化流程[^2]: ```http http://124.71.149.47:5000/smalltools/detail?id=1638074560383299585 ```
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值