草庐IT

【单细胞高级绘图】07.KEGG富集结果展示

TOP生物信息 2023-03-28 原文


这一节画的图是比较新的,图中我用红色箭头标出的是pathway一级注释信息(big annotation,自己想的,非专有名词),纵轴花花绿绿的标注是pathway的二级注释(small annotation)。如何获取注释算一个难点,我上一讲也已经讲过:KEGG通路的从属/注释信息如何获取。

整个图反映的是有多少基因落到了对应的分类里面。

辩证地看,整张图都是pathway注释,没有具体的pathway名称,跟平常做的富集分析很不一样。把图里面的二级注释换成具体的pathway会更好。另外,这个图中横坐标是基因数,但我觉得在富集分析中这个数值并不重要(基因ratio比单个数值重要),我们可以换成p值


在开始画图之前,需要整理一个表格,表格中至少包含:pathway ID、pathway description、pathway注释/big annotation、几个富集指标、被富集到的基因

  • 前三个信息,上一讲的表格已经整理好了,可以拿来用;
  • 几个富集指标clusterprofiler的输出结果中也有;
  • 被富集到的基因貌似clusterprofiler做kegg富集不能直接显示基因symbol,所以返回的结果需要稍微加工一下。

整理表格的代码名称为run.R,如下所示:

library(Seurat)
library(tidyverse)
library(xlsx)

testseu=readRDS("testseu.rds")
Idents(testseu)="anno_new"

### 找差异基因 #########################################################################
marker_celltype=FindAllMarkers(testseu,logfc.threshold = 0.8,only.pos = T)
# 过滤
marker_celltype=marker_celltype%>%filter(p_val_adj < 0.01)
marker_celltype$d=marker_celltype$pct.1-marker_celltype$pct.2
marker_celltype=marker_celltype%>%filter(d > 0.2)
marker_celltype=marker_celltype%>%arrange(cluster,desc(avg_log2FC))
marker_celltype=as.data.frame(marker_celltype)
write.xlsx(marker_celltype,file = "markers_log2fc0.8_padj0.01_pctd0.2.xlsx",row.names = F,col.names = T)

### 富集分析 ###########################################################################
library(clusterProfiler)
library(org.Hs.eg.db)
R.utils::setOption("clusterProfiler.download.method","auto") #https://github.com/YuLab-SMU/clusterProfiler/issues/256

source("syEnrich.R")
syEnrich(marker_celltype,outpath = "markers_log2fc0.8_padj0.01_pctd0.2")

### 挑一类细胞来作为演示 #######################################################
kegg.res=read.xlsx("markers_log2fc0.8_padj0.01_pctd0.2.KEGG.xls",sheetIndex = 1,as.data.frame = T,header = T)
kegg.res=kegg.res%>%filter(p.adjust < 0.05)
kegg.res=kegg.res%>%filter(cluster == "Endothelial")

# 导入上一讲的文件
kegg_info=read.xlsx("kegg_info.xlsx",sheetIndex = 1,startRow = 3)
kegg_info=kegg_info[,c("ID","Pathway","big.annotion")]

# 合并两个表格
kegg.res$ID=str_replace(kegg.res$ID,"hsa","")
kegg.res=kegg.res%>%inner_join(kegg_info,by = "ID")
write.table(kegg.res,file = "kegg.res.txt",quote = F,sep = "\t",row.names = F,col.names = T)

我以单细胞分析中的kegg富集分析作为演示,只取其中一个cluster的富集结果来画图。

上述代码中间用到的富集代码叫syEnrich.R,这个文件只需要输入单细胞seurat对象运行FindAllMarkers得到的差异基因,就可以返回GO/KEGG富集结果,同时被富集到某个通路的基因symbol也会被列出。

运行run.R之后,最终的表格如下图所示:


然后开始画图,代码名称为3.plot.R,这里就不演示了,最终可以得到的图如下:

获取代码

包含这张图会用到的所有代码,数据整理以及画图,超贴心有没有!

这个系列都会采取限时公开的方式共享代码,24小时内是免费的。超过这个时间如何获取,公粽号发送2022A可知。

链接:https://pan.baidu.com/s/1hebbeQH4DgYA8JqFWXO4dg
提取码:abo5

有关【单细胞高级绘图】07.KEGG富集结果展示的更多相关文章

  1. 报告回顾丨模型进化狂飙,DetectGPT能否识别最新模型生成结果? - 2

    导读语言模型给我们的生产生活带来了极大便利,但同时不少人也利用他们从事作弊工作。如何规避这些难辨真伪的文字所产生的负面影响也成为一大难题。在3月9日智源Live第33期活动「DetectGPT:判断文本是否为机器生成的工具」中,主讲人Eric为我们讲解了DetectGPT工作背后的思路——一种基于概率曲率检测的用于检测模型生成文本的工具,它可以帮助我们更好地分辨文章的来源和可信度,对保护信息真实、防止欺诈等方面具有重要意义。本次报告主要围绕其功能,实现和效果等展开。(文末点击“阅读原文”,查看活动回放。)Ericmitchell斯坦福大学计算机系四年级博士生,由ChelseaFinn和Chri

  2. 基于C#实现简易绘图工具【100010177】 - 2

    C#实现简易绘图工具一.引言实验目的:通过制作窗体应用程序(C#画图软件),熟悉基本的窗体设计过程以及控件设计,事件处理等,熟悉使用C#的winform窗体进行绘图的基本步骤,对于面向对象编程有更加深刻的体会.Tutorial任务设计一个具有基本功能的画图软件**·包括简单的新建文件,保存,重新绘图等功能**·实现一些基本图形的绘制,包括铅笔和基本形状等,学习橡皮工具的创建**·设计一个合理舒适的UI界面**注明:你可能需要先了解一些关于winform窗体应用程序绘图的基本知识,以及关于GDI+类和结构的知识二.实验环境Windows系统下的visualstudio2017C#窗体应用程序三.

  3. ruby - gem 推送结果为 "package metadata is missing" - 2

    我正在尝试将我更新的gem推送到ruby​​gems.com并得到以下结果。~/dev/V2/V2GPTI(master)$gembuildv2gpti.gemspecSuccessfullybuiltRubyGemName:v2gptiVersion:0.2File:v2gpti-0.2-universal-darwin-13.gem~/dev/V2/V2GPTI(master)$gempushv2gpti.gemspecERROR:Whileexecutinggem...(Gem::Package::FormatError)packagemetadataismissinginv2g

  4. ruby - 高级语言是否使用数据结构? - 2

    我目前还在上学,正在上一门关于用C++实现数据结构的类(class)。在业余时间,我喜欢使用“高级”语言(主要是Ruby和一些c#)进行编程。既然这些高级语言为你管理内存,你会用数据结构做什么?我可以理解对队列和堆栈的需求,但是您需要在Ruby中使用二叉树吗?还是2-3-4树?为什么?谢谢。 最佳答案 Sosincethesehigherlevellanguagesmanagethememoryforyou,whatwouldyouusedatastructuresfor?使用数据结构的主要原因与垃圾收集无关。但它是以某种方式有效的

  5. ruby - 猴子修补 float 中缀运算符产生意想不到的结果 - 2

    重新定义Float#/似乎没有效果:classFloatdef/(other)"magic!"endendputs10.0/2.0#=>5.0但是当另一个中缀运算符Float#*被重新定义时,Float#/突然采用了新的定义:classFloatdef/(other)"magic!"enddef*(other)"spooky"endendputs10.0/2.0#=>"magic!"我很想知道是否有人可以解释这种行为的来源,以及其他人是否得到相同的结果。ruby:ruby2.0.0p353(2013-11-22)[x64-mingw32]要快速确认错误,请运行thisscript.

  6. ruby-on-rails - 尝试编辑时,Rails form_for 结果为 POST 而不是 PUT - 2

    我正在使用Rails4并遇到以下错误。RoutingErrorNoroutematches[POST]"/logs/1/meals/13/edit我正在使用:meal传递模型对象的form_for,并且编辑页面正确呈现。但是,Rails似乎并没有检查膳食对象是否已经保存,因此它一直尝试将表单发送到#create操作并尝试发出POST请求,而不是将表单发送到更新操作并进行当我点击提交时一个PUT请求。我如何让form_for识别我正在尝试更新现有对象并且需要PUT而不是POST?其他一切正常,我已经运行了所有迁移。我是Rails的新手,几乎一整天都在尝试自己解决这个问题。请帮忙!请注意,

  7. ruby - 为什么会.is_a?和 .class 给出相互矛盾的结果? - 2

    我有三个属于同一个类的对象。一个是通过Item.new创建的,另外两个是从数据库(Mongoid)中提取的。我将这些对象中的一个/任何一个传递给另一个方法,并通过is_a?检查该方法中的类型:definitialize(item,attrs=nil,options=nil)super(attrs,options)raise'invaliditemobject'unlessitem.is_a?(Item)好吧,这次加薪被击中了。所以我在Rails控制台中检查类、is_a和instance_of。我得到相互矛盾的结果。为什么它们有相同的class但只有其中一个是那个class的instan

  8. ruby - 带有位串的意外解包结果 - 2

    为什么当我打开irb并运行时放'A'.unpack("B8")我得到01000001但是当我运行放'A'.unpack("B4B4")我只得到0100而不是[0100,0001]?unpack的分辨率是不是只有一个完整的字节?一点都不差? 最佳答案 让我们做一些测试来理解行为:>'A'.unpack('B8')=>["01000001"]它返回char'A'的8个最高有效位(MSB)>'A'.unpack('B4')=>["0100"]它返回char'A'的4MSBs>'A'.unpack('B16')=>["01000001"]它

  9. ruby - 为什么我会看到这两个几乎相同的 Ruby 正则表达式模式的不同结果,为什么一个匹配我认为不应该匹配的内容? - 2

    使用Ruby1.9.2,我在IRB中有以下Ruby代码:>r1=/^(?=.*[\d])(?=.*[\W]).{8,20}$/i>r2=/^(?=.*\d)(?=.*\W).{8,20}$/i>a=["password","1password","password1","pass1word","password1"]>a.each{|p|puts"r1:#{r1.match(p)?"+":"-"}\"#{p}\"".ljust(25)+"r2:#{r2.match(p)?"+":"-"}\"#{p}\""}这会产生以下输出:r1:-"password"r2:-"password"r1:

  10. ruby - 为什么 Gemfile 语义版本控制运算符 (~>) 会产生与一个数字不一致的结果? - 2

    gemspec语义版本控制运算符~>(又名twiddle-wakka,又名pessimistic运算符)允许限制gem版本但允许进行一些升级。我经常看到它可以读作:"~>3.1"=>"Anyversion3.x,butatleast3.1""~>3.1.1"=>"Anyversion3.1.x,butatleast3.1.1"但是有了一个数字,这条规则就失效了:"~>3"=>"Anyversionx,butatleast3"*NOTTRUE!*"~>3"=>"Anyversion3.x"*True.Butwhy?*如果我想要“任何版本3.x”,我可以只使用“~>3.0”,这是一致的。就

随机推荐