基础生信分析——富集分析(1)

本文介绍了如何在基因差异表达分析后使用KEGG和GO数据库进行富集分析,包括富集分析的概念、步骤,以及如何在R中使用`clusterProfiler`包进行环境设置、数据读取、统计分析和可视化结果。重点讲解了如何进行环境配置、导入DEGs数据、进行GO和KEGG数据库的分析,并展示了基础的柱状图和气泡图可视化方法。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

基础生信分析——富集分析(1)



前言

KEGG(Kyoto Encyclopedia of Genes and Genomes)和GO(Gene Ontology)是两种常用的基因功能和通路富集分析方法,可以帮助研究人员理解基因集合中的功能和通路特征。以下是对KEGG/GO富集分析的简要介绍:

KEGG富集分析:KEGG富集分析基于KEGG数据库,对基因集合中的基因进行功能注释和通路富集分析。通过将基因集合与KEGG数据库中的通路注释进行比较,可以确定哪些通路在给定基因集合中过度表示。这有助于理解基因集合中的生物学过程、信号传导通路或代谢通路的重要性。

GO富集分析:GO富集分析基于Gene Ontology,它是一个标准化的基因和蛋白质功能注释系统。GO由三个纬度组成:分子功能(Molecular Function)、细胞组分(Cellular Component)和生物过程(Biological Process)。GO富集分析可以帮助确定给定基因集合中具有显著富集的功能或过程。

这两种富集分析方法通常包括以下步骤:

基因集定义:首先,需要明确你感兴趣的基因集合,可以是差异表达基因、某个生物学过程相关的基因等。

统计显著性分析:在富集分析中,使用统计方法(如超几何分布、Fisher’s精确检验等)来衡量基因集与注释数据库中的功能或通路之间的关联性,并确定哪些功能或通路在基因集中过度表示。

多重假设检验校正:由于分析过程中进行了多次比较,通常需要进行多重假设检验校正,例如Bonferroni、Benjamini-Hochberg调整等。

结果解释:根据分析结果,可以识别出在给定基因集合中显著富集的KEGG通路或GO功能,并进一步解释这些富集结果的生物学意义。

KEGG/GO富集分析是生物学研究中常用的功能注释和通路分析方法,可以帮助研究人员洞察基因集背后的生物学特性和功能,从而深入理解基因集的生物学意义。


一、富集分析是什么?

简单来说,差异分析之后我们仍然获得非常多的DEGs,我们需要进一步对这些DEGs做一些功能注释或者对某些特点功能和通路的基因做深度筛选,这时候就需要进行富集分析,也叫功能分析

二、使用步骤

1.环境

# 下载所需要的R包 ----------------------------------------------------------------
if (!require("BiocManager", quietly = TRUE))
  install.packages("BiocManager")
# 检查 clusterProfiler 包是否已经安装
if (!require("clusterProfiler", quietly = TRUE)) {
  # 如果未安装,则使用 BiocManager 安装
  if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
  BiocManager::install("clusterProfiler")
}
# 检查 org.Hs.eg.db 包是否已经安装
if (!require("org.Hs.eg.db", quietly = TRUE)) {
  # 如果未安装,则使用 BiocManager 安装
  if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
  BiocManager::install("org.Hs.eg.db")
}
# 检查 ggplot2 包是否已经安装
if (!require("ggplot2", quietly = TRUE)) {
  # 如果未安装,则使用 install.packages 函数安装
  install.packages("ggplot2")
}
# 检查 dplyr 包是否已经安装
if (!require("dplyr", quietly = TRUE)) {
  # 如果未安装,则使用 install.packages 函数安装
  install.packages("dplyr")
}
# 检查 ggsci 包是否已经安装
if (!require("ggsci", quietly = TRUE)) {
  # 如果未安装,则使用 install.packages 函数安装
  install.packages("ggsci")
}
# 检查 GOplot 包是否已经安装
if (!require("GOplot", quietly = TRUE)) {
  # 如果未安装,则使用 install.packages 函数安装
  install.packages("GOplot")
}
# 加载所需要的R包 ----------------------------------------------------------------
library(clusterProfiler) #富集分析主要的包
library(org.Hs.eg.db)#查找物种注释信息
library(ggplot2)#分面绘图所需
library(dplyr)#整理数据
library(ggsci)#配色
library(GOplot)#可视化

2.读入DEGs数据

# 读取数据:基因差异分析的结果 ----------------------------------------------------------
Gene <- read.csv("DEG_deseq2.csv",row.names = 1)
Genes <- bitr(Gene$Gene_symbol,  
              fromType = "SYMBOL", #输入数据的类型
              toType = c("ENTREZID"), #要转换的数据类型
              OrgDb = org.Hs.eg.db) #物种

3.GO数据库

分析计算

# GO富集分析 --------------------------------------------------------------------
GO <- enrichGO(gene = Genes$ENTREZID, #输入基因的"ENTREZID"
               OrgDb = org.Hs.eg.db,#注释信息
               keyType = "ENTREZID",
               ont = "all",     #可选条目BP/CC/MF
               pAdjustMethod = "BH", #p值的校正方式
               pvalueCutoff = 1,   #pvalue的阈值
               qvalueCutoff = 1, #qvalue的阈值
               minGSSize = 5,
               maxGSSize = 5000,
               readable = TRUE)   #是否将entrez id转换为symbol
GO_result <- as.data.frame(GO)#转化结果

可视化

# 基础柱状图 -------------------------------------------------------------------
barplot(GO)#单个条目
barplot(GO, drop = TRUE, showCategory =6,split="ONTOLOGY") +
facet_grid(ONTOLOGY~.,scale='free')#多个条目

# 美化版本 --------------------------------------------------------------------
#挑选每个分类中前10显著的条目
BP_top10 <- GO_result %>%
  filter(ONTOLOGY == "BP") %>%  # 筛选 ONTOLOGY 为 BP 的子集
  arrange(p.adjust) %>%        # 按照 p.adjust 的升序排列
  head(10)                      # 提取前 10行

CC_top10 <- GO_result %>%
  filter(ONTOLOGY == "CC") %>%  # 筛选 ONTOLOGY 为 CC 的子集
  arrange(p.adjust) %>%        # 按照 p.adjust 的升序排列
  head(10)                      # 提取前 10行

MF_top10 <- GO_result %>%
  filter(ONTOLOGY == "MF") %>%  # 筛选 ONTOLOGY 为 MF 的子集
  arrange(p.adjust) %>%        # 按照 p.adjust 的升序排列
  head(10)                      # 提取前 10行

merge_data <- rbind(BP_top10,CC_top10,MF_top10)#合并三个条目
merge_data$ONTOLOGY <- factor(merge_data$ONTOLOGY, levels = c("BP", "CC","MF"))
merge_data$logPvalue <- (-log(merge_data$p.adjust))#数据处理

colnames(merge_data)#检查列名
merge_data$Description <- factor(merge_data$Description, levels = unique(merge_data$Description[order(merge_data$ONTOLOGY)]))

ggplot(merge_data, aes(x = Description, y = logPvalue, fill = ONTOLOGY)) + #指定数据为merge_data,以Description作为x轴,以logPvalue作为y轴,以ONTOLOGY作为颜色填充。
  geom_bar(stat = "identity", width = 0.8, color = "black") + #添加一个柱状图层,使用identity统计方法(不进行任何变换),每个柱子的宽度为0.8,颜色为黑色
  scale_x_discrete(limits = unique(merge_data$Description[order(merge_data$ONTOLOGY)])) + #设置x轴的标签顺序,根据ONTOLOGY对Description排序,以便按照ONTOLOGY的顺序显示。unique()函数保证每个标签只出现一次。
  coord_flip() + #翻转坐标轴,交换x轴和y轴,使得柱状图变为横向。
  scale_y_continuous(expand = c(0,0))+ #设置y轴的范围为连续型变量(continuous),设置expand参数为c(0,0),意味着y轴不需要扩展(不需要增加额外的空白)。
  theme(panel.background = element_rect(fill = "white"))+ #设置背景为白色,使用element_rect()函数设置元素的矩形形状
  theme(axis.text = element_text(size = 16, family = "Arial"),
        axis.title = element_text(size=18,family="Arial",colour = "black"),
        legend.text = element_text(size = 16, family = "Arial"),
        legend.title = element_text(size = 18, family = "Arial"),
        legend.position = "top")+ #设置各种文本元素的字体、大小和颜色,包括坐标轴文本、坐标轴标题、图例文本和图例标题。将图例放置在顶部。
  labs(x = "GO term", y = "-log P-value")+
  scale_fill_nejm(alpha = 0.8) #使用NEJM(New England Journal of Medicine)颜色调色板,设置颜色透明度为0.8,对填充颜色进行比例缩放。

在这里插入图片描述

3.KEGG数据库

分析计算

# KEGG富集分析----------------------------------------------------------------
KEGG <- enrichKEGG(gene = Genes$ENTREZID,
                   organism = "hsa", 
                   keyType = "kegg", #KEGG数据库
                   pAdjustMethod = "BH",
                   pvalueCutoff = 1,
                   qvalueCutoff = 1)

可视化

#柱状图----------------------------------------------------------------------
barplot(KEGG,
       x = "GeneRatio",
       color = "p.adjust",
       showCategory = 10,
       title = "KEGG_enrichment")
# 气泡图----------------------------------------------------------------------
#根据基因count数排序
dotplot(KEGG)
#根据P值排序
df <- as.data.frame(KEGG)#转换数据
df_sorted <- df[order(df$p.adjust), ][1:10, ]#前10名
colnames(df_sorted)
df_sorted$Description
df_sorted$Description <- gsub(" - Homo sapiens \\(human\\)", "", df_sorted$Description)

ggplot(df_sorted, aes(x=GeneRatio, y=Description, size=-log10(p.adjust), color=p.adjust)) + ## 生成散点图,以GeneRatio为x轴,Description为y轴,-log10(p.adjust)为点的大小,p.adjust为点的颜色
 geom_point(alpha=0.7) + # 绘制散点图层,设置透明度为0.7
 scale_color_gradient(low="blue", high="red") +
 labs(title="The most enrichment KEGG pathway", x="GeneRatio", y="KEGG pathway", size="-log10(P-value)", color="P-value") + #设置颜色映射为渐变色,低值为蓝色,高值为红色
 theme_minimal() + #设置图表的标题、x轴和y轴标签,以及点的大小和颜色的图例标签
 theme(axis.text.y = element_text(face = 'bold',size = 12),
       axis.title = element_text(size = 15),
       plot.title = element_text(size = 18, face = "bold")) #设置主题为最小化主题,调整y轴的文本大小和字体粗细,以及标题、轴标签的字体大小

在这里插入图片描述


总结

以上是简单的对DEGs做GO/KEGG富集分析

### 如何使用 PlantTFDB 下载基因家族数据并进行息学分析 #### 获取所需工具和资源 为了有效地从 PlantTFDB 数据库获取基因家族的数据,并对其进行后续的息学分析,需要准备一些基本工具和软件环境。这通常包括但不限于 Python 或 R 编程语言及其相应的包来处理序列文件以及执行统计测试。 #### 访问 PlantTFDB 并定位目标基因家族 访问 [PlantTFDB](http://planttfdb.cbi.pku.edu.cn/) 官方网站,在首页可以找到不同版本之间的切换选项。选择最新发布的版本以获得最全面的息。通过导航栏中的“Browse TFs by Family”,可以选择特定植物物种下的各个转录因子家族列表[^1]。 #### 下载感兴趣的基因家族成员 一旦确定了要研究的具体基因家族,点击进入该页面后可以看到详细的成员清单。这里提供了多种格式供下载,比如 FASTA 序列、GFF 文件等。对于大多数情况来说,FASTA 格式的核苷酸或氨基酸序列是最常用的输入形式之一用于下游的各种计算物学流程。 #### 初步探索与预处理所获资料 利用命令行工具如 `grep` 和 `awk` 可快速筛选出关注区域内的记录;而像 Biopython 这样的库则有助于解析复杂的多序列比对结果或是构建进化树模型。下面给出一段简单的Python脚本作为例子展示如何读取本地保存下来的fasta文档: ```python from Bio import SeqIO def read_fasta(file_path): sequences = [] with open(file_path, 'r') as handle: for record in SeqIO.parse(handle, "fasta"): seq_id = record.id sequence = str(record.seq) sequences.append((seq_id, sequence)) return sequences # 使用函数加载 fasta 文件 gene_family_sequences = read_fasta('path_to_your_downloaded_file.fasta') print(f'Total number of genes loaded: {len(gene_family_sequences)}') ``` 此段代码实现了从指定路径下读入 FASTA 文件的功能,并将其转换成易于操作的形式——即由元组组成的列表,每个元组包含两个元素分别是基因ID与其对应的核酸/蛋白序列字符串表示法。 #### 执行功能富集分析和其他高级应用 当拥有了足够的原始材料之后就可以着手开展更深入的研究工作了。例如,可以通过在线平台 DAVID 来实施 GO term 富集检验从而揭示潜在的作用机制;或者借助 STRING 数据库查询已知互作伙伴进而推测未知关联关系等等[^5]。
评论 3
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

柳叶刀和鼠标

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

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

余额充值