草庐IT

转录组分析(三)--GO及KEGG富集分析及可视化(包括条形图及气泡图)(对差异基因的可视化)

千万别加香菜 2023-03-28 原文
#### 安装包
install.packages("BiocManager")
BiocManager::install(version = "3.13")
BiocManager::install("gprofiler2")
BiocManager::install("clusterProfiler")
BiocManager::install("AnnotationHub")
BiocManager::install("org.Bt.eg.db")

GO分析(上下调基因分开做,可用于BP,CC,MF分开画图)

## 方法2:下载到本地加载,每次使用上传,(推荐)
library(AnnotationDbi)
setwd("D:/生信/sheep.db")
sheep.db <- AnnotationDbi::loadDb('org.Ovis_aries.eg.sqlite')
columns(sheep.db)

## 牛
library(AnnotationDbi)
setwd("D:/生信/cattle.db")
bos.db <- AnnotationDbi::loadDb('org.Bt.eg.db.sqlite')
columns(bos.db)
#读取基因列表文件中的基因名称
genes <- read.table('GENEID.txt', header = T, stringsAsFactors = FALSE,sep = '\t')
library("org.Bt.eg.db")
library(clusterProfiler)
#GO富集分析,上下调分开做
#对于加载的注释库的使用,以上述为例,就直接在 OrgDb 中指定人(org.Hs.eg.db)或绵羊(sheep)
enrich.go <- enrichGO(gene = genes$V1,  #基因列表文件中的基因名称
                      OrgDb = "org.Bt.eg.db",  #指定物种的基因数据库
                      keyType = 'SYMBOL',  #指定给定的基因名称类型,例如这里以 SYMBOL 为例
                      ont = 'CC',  #可选 BP、MF、CC,也可以指定 ALL 同时计算 3 者
                      pAdjustMethod = 'BH',  #指定 p 值校正方法
                      pvalueCutoff = 1,  #指定 p 值阈值,不显著的值将不显示在结果中
                      qvalueCutoff = 1  #指定 q 值阈值,不显著的值将不显示在结果中
                     )
# 例如上述指定 ALL 同时计算 BP、MF、CC,这里将它们作个拆分后输出
BP <- enrich.go[enrich.go$ONTOLOGY=='BP', ]
CC <- enrich.go[enrich.go$ONTOLOGY=='CC', ]
MF <- enrich.go[enrich.go$ONTOLOGY=='MF', ]
write.table(as.data.frame(BP), 'HUgo.BP.txt', sep = '\t', row.names = FALSE, quote = FALSE)
write.table(as.data.frame(CC), 'HUgo.CC.txt', sep = '\t', row.names = FALSE, quote = FALSE)
write.table(as.data.frame(MF), 'HUgo.MF.txt', sep = '\t', row.names = FALSE, quote = FALSE)

GO分析(上下调一起做,强烈推荐,可以看整体结果,并可选取其中的部分画图)

#up及down基因的合并处理 
gene_down <- read.table('up_gene.txt', header = T, stringsAsFactors = FALSE,sep = '\t')
gene_up <- read.table('down_gene.txt', header = T, stringsAsFactors = FALSE,sep = '\t')
# 转基因ID
library(clusterProfiler)
deg.up.gene <- bitr(gene_up$x, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = bos.db)$ENTREZID 
deg.down.gene <- bitr(gene_down$x, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = bos.db)$ENTREZID
# 将up 和 down基因合并到一个表格
deg.ei.ls <- list(up = deg.up.gene, down = deg.down.gene)

# GO
library("org.Bt.eg.db")
GO.cmp <- compareCluster(deg.ei.ls,
                         fun = "enrichGO",
                         qvalueCutoff = 0.1, 
                         pvalueCutoff = 0.01,
                         pAdjustMethod = 'BH',
                         ont = 'BP', #可选 BP、MF、CC,也可以指定 ALL 同时计算 3 者,建议分开做输出,默认按照pvalue升序排列的
                         OrgDb = "org.Bt.eg.db")
write.table(GO.cmp@compareClusterResult,'GO-up-down.txt',sep = '\t',quote = F)

GO 条形图(上下调一起,不过BP,CC,MF等需要分开做)

library(ggplot2)
library(ggsci)
ddf<-read.table('GO-up-down.txt',header = T,sep = '\t')
dat=ddf
dat$Cluster = ifelse(dat$Cluster =='up',1,-1)
dat$pvalue = -log10(dat$pvalue)
dat$pvalue <- dat$pvalue*dat$Cluster
dat=dat[order(dat$pvalue,decreasing = F),]

ggplot(dat, aes(x=reorder(Description,order(pvalue, decreasing = F)), 
                y=pvalue,fill=Cluster)) + 
  geom_bar(stat="identity",width = 0.78,position = position_dodge(0.7)) + 
  scale_fill_gradient(low="#3685af",high="#af4236",guide = FALSE)+
  scale_x_discrete(name ="GO Terms") +
  #scale_y_discrete(labels =Description )
  scale_y_continuous(name ="-log10Pvalue",expand = c(0,0)) +
  coord_flip() +  # 翻转坐标
  theme(panel.grid.major.y = element_blank(),panel.grid.minor.y = element_blank())+ #删除纵向网络
  theme(axis.text.y = element_blank())+ #刻度标签空白
  theme(axis.ticks = element_blank())+ #刻度为空白
  theme(panel.border = element_blank())+ #边框为空白
  theme(axis.text=element_text(face = "bold",size = 15),
        axis.title = element_text(face = 'bold',size = 15), #坐标轴字体
        plot.title = element_text(size = 20,hjust = 0.3), 
        legend.position = "top",panel.grid = element_line(colour = 'white'))+
  geom_text(aes(label=Description),size=3.5,hjust=ifelse(dat$pvalue>0,1.5,-0.5))+ # 两侧加上标签文字
  ggtitle("The most enriched GO-BP terms")  # 添加标题

GO气泡图(可上下调一起做,也不用将BP,CC,MF分开,最好选择适当数量)

从这位老哥那里薅过来的,这里是原文链接(https://blog.csdn.net/sweet_yemen/article/details/125457274)

# 富集气泡图
 # BiocManager::install("openxlsx")
library(ggplot2)#这个是一个非常厉害的绘图R包
library(openxlsx)#这个包是用来导入格式为XLSX的数据的

goinput <- read.xlsx("go-barplot.xlsx")#载入GO富集分析的结果文件
# 这里强调一下富集后的GeneRatio是以n/n这种形式的,需要换成小数形式,在excel中使用=--("0 "&A1)
# A1代表对应的单元格
# 设置XY轴
x=goinput$GeneRatio#设置以哪个数据为X轴
y=factor(goinput$Description,levels = goinput$Description)#设置Y轴
# 绘制画布
p = ggplot(goinput,aes(x,y))#开始以X轴和Y轴绘制画布
p#这里可以查看一下绘制的程度了,学代码时多看看每个新的对象有助于理解代码

#绘图,大概的意思是,又设置了除XY轴数据(Ratio,GOterm)之后又设置了三个维度的数据。
#三个数据分别是以count为点的大小,以颜色的由嫩绿到深粉红代表P值的对数值,以点的性状代表GO类型(MF,CC,BP)
p1 = p + geom_point(aes(size=Count,color=-0.5*log(pvalue),shape=ONTOLOGY,))+
  scale_color_gradient(low = "SpringGreen", high = "DeepPink")
p1
#添加各类标签
p2 = p1 + labs(color=expression(-log[10](pvalue)),
               size="Count",
               x="Ratio",
               y="Go_term",
               title="Go enrichment of test Genes")
p2

KEGG, 上下调一起做(适合看上下调整体的结果)

#up及down基因的合并处理 
gene_down <- read.table('up_gene.txt', header = T, stringsAsFactors = FALSE,sep = '\t')
gene_up <- read.table('down_gene.txt', header = T, stringsAsFactors = FALSE,sep = '\t')
# 转基因ID
library(clusterProfiler)
deg.up.gene <- bitr(gene_up$x, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = bos.db)$ENTREZID 
deg.down.gene <- bitr(gene_down$x, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = bos.db)$ENTREZID
# 将up 和 down基因合并到一个表格
deg.ei.ls <- list(up = deg.up.gene, down = deg.down.gene)

library(clusterProfiler)
R.utils::setOption("clusterProfiler.download.method",'auto') ##一定要用这个,不然报错别怪我哦
kegg.cmp <- compareCluster(deg.ei.ls,
                           fun = "enrichKEGG",
                           qvalueCutoff = 1, 
                           pvalueCutoff = 1,
                           pAdjustMethod = 'BH',
                           organism = "bta")
dotplot(kegg.cmp)
write.table(kegg.cmp@compareClusterResult,'kegg-up-down.txt',sep = '\t',quote = F)

kegg, 上下调分开做(适合画图)

library(clusterProfiler)
genes <- read.table('txt', header = T, stringsAsFactors = FALSE,sep = '\t')
covID <- bitr(genes$V1, fromType = "SYMBOL",
                   toType="ENTREZID",OrgDb= bos.db , drop = TRUE) # 转换为entrezID,KEGG只识别这种ID,报错请检查自己输入的文件,对自己的代码有自信
write.table(covID,'ENTREZID-up.txt',sep = '\t',quote = F) # 写出去保存,也可以不保存,看个人需求

#每次打开R计算时,它会自动连接kegg官网获得最近的物种注释信息,因此数据库一定都是最新的
search_kegg_organism("taurus", by="scientific_name") # 寻找对应动物的organism

R.utils::setOption("clusterProfiler.download.method",'auto') ##一定要用这个,不然报错别怪我
enrich_kegg <- enrichKEGG(gene = covID$ENTREZID, # 基因列表文件中的基因名
                          organism = 'bta',  # 指定物种,牛为bta
                          keyType = 'kegg', # kegg富集
                          pAdjustMethod = 'BH', # 指定矫正方法
                          pvalueCutoff = 1, # 指定p值阈值,不显著的值将不显示在结果中
                          qvalueCutoff = 1) # 指定q值阈值,不显著的值将不显示在结果中
#输出结果
write.table(enrich_kegg@result,'kegg-down.txt',quote = F,sep = '\t')

KEGG富集气泡图(上下调分开做推荐)

library(ggplot2)
dotplot(enrich_kegg,
        showCategory = 15,#展示前15个点,默认为10个
        color = 'pvalue',
        title = "kegg_dotplot")

KEGG富集条形图

barplot(enrich_kegg,       
        showCategory = 15,#展示前15个点,默认为10个
        color = 'pvalue',
        title = "kegg_dotplot")

有关转录组分析(三)--GO及KEGG富集分析及可视化(包括条形图及气泡图)(对差异基因的可视化)的更多相关文章

  1. ruby - Ruby 中的波形可视化 - 2

    我即将开始一个将录制和编辑音频文件的项目,我正在寻找一个好的库(最好是Ruby,但会考虑Java或.NET以外的任何库)以进行实时可视化波形。有人知道我应该从哪里开始搜索吗? 最佳答案 要流入浏览器的数据量很大。Flash或Flex图表可能是唯一能提高内存效率的解决方案。Javascript图表往往会因大型数据集而崩溃。 关于ruby-Ruby中的波形可视化,我们在StackOverflow上找到一个类似的问题: https://stackoverflow.c

  2. ruby - 在 StockChart (highchart) 中以编程方式显示柱形图的工具提示 - 2

    我有一个Highstock图表(带有标记和阴影的线条),并且想以编程方式显示一个highstock工具提示,例如,当我选择某个表上的一行(包含图表数据)我想显示相应的highstock工具提示。这可能吗? 最佳答案 股票图表thissolution不起作用:在thisexample你必须更换这个:chart.tooltip.refresh(chart.series[0].data[i]);为此:chart.tooltip.refresh([chart.series[0].points[i]]);解决方案可用here.

  3. 建模分析 | 平面2R机器人(二连杆)运动学与动力学建模(附Matlab仿真) - 2

    目录0专栏介绍1平面2R机器人概述2运动学建模2.1正运动学模型2.2逆运动学模型2.3机器人运动学仿真3动力学建模3.1计算动能3.2势能计算与动力学方程3.3动力学仿真0专栏介绍?附C++/Python/Matlab全套代码?课程设计、毕业设计、创新竞赛必备!详细介绍全局规划(图搜索、采样法、智能算法等);局部规划(DWA、APF等);曲线优化(贝塞尔曲线、B样条曲线等)。?详情:图解自动驾驶中的运动规划(MotionPlanning),附几十种规划算法1平面2R机器人概述如图1所示为本文的研究本体——平面2R机器人。对参数进行如下定义:机器人广义坐标

  4. 网站日志分析软件--让网站日志分析工作变得更简单 - 2

    网站的日志分析,是seo优化不可忽视的一门功课,但网站越大,每天产生的日志就越大,大站一天都可以产生几个G的网站日志,如果光靠肉眼去分析,那可能看到猴年马月都看不完,因此借助网站日志分析工具去分析网站日志,那将会使网站日志分析工作变得更简单。下面推荐两款网站日志分析软件。第一款:逆火网站日志分析器逆火网站日志分析器是一款功能全面的网站服务器日志分析软件。通过分析网站的日志文件,不仅能够精准的知道网站的访问量、网站的访问来源,网站的广告点击,访客的地区统计,搜索引擎关键字查询等,还能够一次性分析多个网站的日志文件,让你轻松管理网站。逆火网站日志分析器下载地址:https://pan.baidu.

  5. ABB-IRB-1200运动学分析MATLAB RVC工具分析+Simulink-Adams联合仿真 - 2

    一、机器人介绍        此处是基于MATLABRVC工具箱,对ABB-IRB-1200型号的微型机械臂进行正逆向运动学分析,并利Simulink工具实现对机械臂进行具有动力学参数的末端轨迹规划仿真,最后根据机械模型设计Simulink-Adams联合仿真。 图1.ABBIRB 1200尺寸参数示意图ABBIRB 1200提供的两种型号广泛适用于各作业,且两者间零部件通用,两种型号的工作范围分别为700 mm 和 900 mm,大有效负载分别为 7 kg 和5 kg。 IRB 1200 能够在狭小空间内能发挥其工作范围与性能优势,具有全新的设计、小型化的体积、高效的性能、易于集成、便捷的接

  6. 关于Qt程序打包后运行库依赖的常见问题分析及解决方法 - 2

    目录一.大致如下常见问题:(1)找不到程序所依赖的Qt库version`Qt_5'notfound(requiredby(2)CouldnotLoadtheQtplatformplugin"xcb"in""eventhoughitwasfound(3)打包到在不同的linux系统下,或者打包到高版本的相同系统下,运行程序时,直接提示段错误即segmentationfault,或者Illegalinstruction(coredumped)非法指令(4)ldd应用程序或者库,查看运行所依赖的库时,直接报段错误二.问题逐个分析,得出解决方法:(1)找不到程序所依赖的Qt库version`Qt_5'

  7. ruby-on-rails - 如何使用 ruby​​-prof 和 JMeter 分析 Rails - 2

    我想使用ruby​​-prof和JMeter分析Rails应用程序。我对分析特定Controller/操作/或模型方法的建议方法不感兴趣,我想分析完整堆栈,从上到下。所以我运行这样的东西:RAILS_ENV=productionruby-prof-fprof.outscript/server>/dev/null然后我在上面运行我的JMeter测试计划。然而,问题是使用CTRL+C或SIGKILL中断它也会在ruby​​-prof可以写入任何输出之前杀死它。如何在不中断ruby​​-prof的情况下停止mongrel服务器? 最佳答案

  8. Unity数据可视化图表插件XCharts3.0发布 - 2

    Unity数据可视化图表插件XCharts3.0发布历时8个多月,业余时间,断断续续,XCharts3.0总算发布了。如果要打个满意度,我给3.0版本来个80分。对于代码框架结构设计的调整改动,基本符合预期,甚是满意。相比之前的1.0和2.0版本,我认为3.0才是一个拿得出手给广大开发者使用的版本。1.0发布的时候,很兴奋,从0.1到1.0,也磨了一年,真的等不及想给大家试用了,还特地写过一篇文章以示庆祝。那个时候,1.0虽然还还不够完善,功能也不够丰富,但它是XCharts的开始,没有1.0,也就没有后面的2.0和3.0。后面的2.0发布,做了很多改进和优化,随着版本迭代,慢慢的发现有不少硬

  9. 【Unity游戏破解】外挂原理分析 - 2

    文章目录认识unity打包目录结构游戏逆向流程Unity游戏攻击面可被攻击原因mono的打包建议方案锁血飞天无限金币攻击力翻倍以上统称内存挂透视自瞄压枪瞬移内购破解Unity游戏防御开发时注意数据安全接入第三方反作弊系统外挂检测思路狠人自爆实战查看目录结构用il2cppdumper例子2-森林whoishe后记认识unity打包目录结构dll一般很大,因为里面是所有的游戏功能编译成的二进制码游戏逆向流程开发人员代码被编译打包到GameAssembly.dll中使用il2ppDumper工具,并借助游戏名_Data\il2cpp_data\Metadata\global-metadata.dat

  10. 驱动开发:内核无痕隐藏自身分析 - 2

    在笔者前面有一篇文章《驱动开发:断链隐藏驱动程序自身》通过摘除驱动的链表实现了断链隐藏自身的目的,但此方法恢复时会触发PG会蓝屏,偶然间在网上找到了一个作者介绍的一种方法,觉得有必要详细分析一下他是如何实现的进程隐藏的,总体来说作者的思路是最终寻找到MiProcessLoaderEntry的入口地址,该函数的作用是将驱动信息加入链表和移除链表,运用这个函数即可动态处理驱动的添加和移除问题。MiProcessLoaderEntry(pDriverObject->DriverSection,1)添加MiProcessLoaderEntry(pDriverObject->DriverSection,

随机推荐