自定义背景基因非模式生物Go富集：使用uniprot数据自行建库，不需要orgdb
需要做一个非模式生物的Go富集，在目前能搜索到的教程里面，大多是通过AnnotationHub检索并抓取OrgDb。但是本次发现这个方法获得的数据库有点问题，
具体过程是下载的org.db和扒下来的某物种全基因组交集太小，导致最后go富集的时候p值太大
按说go富集本来p值就是很容易小的，结果做出的结果0.01都费劲
于是本次发现可以通过uniprot数据库如今强大的注释信息来跳过复杂的过程，直接搞定go富集
这种方法不止可以应用于某物种的独立研究，可以任意设计自己的背景基因信息
首先回顾一下org.db依赖的go富集方法
1.Org.db依赖的GO富集
以鸭子(Anas,Taxa_id:8835)为例
library(clusterProfiler)# 3.2版本以上的R需要biomanager安装
library(biomaRt)
library(AnnotationHub)
hub <- AnnotationHub(localHub=FALSE)
#AnnotationHub在最开始使用的时候会建立本地的缓存
anas_result<-query(hub,'Anas')
unique(anas_result$rdataclass)>output:
[1] "GRanges" "Inparanoid8Db" "TwoBitFile" "ChainFile" "EnsDb" "SQLiteFile"
[7] "OrgDb"
#这里显示了annotationHub所储存到不同类别的数据库anas_result<-anas_result[anas_result$rdataclass=='OrgDb']#提取OrgDb类型
anas_result$description>output:[1] "NCBI gene ID based annotations about Ananas comosus" [2] "NCBI gene ID based annotations about Ananas comosus_var._comosus" [3] "NCBI gene ID based annotations about Ananas lucidus" [4] "NCBI gene ID based annotations about Anas boschas" [5] "NCBI gene ID based annotations about Anas domesticus" [6] "NCBI gene ID based annotations about Anas platyrhynchos_f._domestica"[7] "NCBI gene ID based annotations about Anas platyrhynchos" [8] "NCBI gene ID based annotations about Anas olor" [9] "NCBI gene ID based annotations about Anas atrata"
[10] "NCBI gene ID based annotations about Anas fuligula"
[11] "NCBI gene ID based annotations about Drosophila ananassae"
[12] "NCBI gene ID based annotations about Drosophila annanassae"
#这里就是显示了所收录的orgdb，其中6到10都是鸭子，具体可以调取其中的内容来看分类号，产生的时间anas_result[6]>output:
'''
AnnotationHub with 1 record
# snapshotDate(): 2022-04-25
# names(): AH101743
# $dataprovider: ftp://ftp.ncbi.nlm.nih.gov/gene/DATA/
# $species: Anas platyrhynchos_f._domestica
# $rdataclass: OrgDb
# $rdatadateadded: 2022-04-21
# $title: org.Anas_platyrhynchos_f._domestica.eg.sqlite
# $description: NCBI gene ID based annotations about Anas platyrhynchos_f._domestica
# $taxonomyid: 8839
# $genome: NCBI genomes
# $sourcetype: NCBI/UniProt
# $sourceurl: ftp://ftp.ncbi.nlm.nih.gov/gene/DATA/, ftp://ftp.uniprot.org/pub/databases/uniprot/current_rele...
# $sourcesize: NA
# $tags: c("NCBI", "Gene", "Annotation")
# retrieve record with 'object[["AH101743"]]'
'''到了这一步，我们算是找到了可以用的数据库，之后用这个来做go富集就OK
anas_org<-anas_result[["AH101742"]]
随后开始进行富集分析，gene_symbol是已经导入进去的基因名
gene_symbol = mapIds(x = anas_org,keys = gene_symbol,keytype = "ENTREZID",column = 'PMID')erich.go.BP = enrichGO(gene = GO,OrgDb = anas_org,keyType ="GO",ont = "BP",pvalueCutoff =0.05, qvalueCutoff = 0.05) erich.go.CC = enrichGO(gene = GO,OrgDb = anas_org,keyType ="GO",ont = "CC",pvalueCutoff = 0.05,qvalueCutoff = 0.05)erich.go.MF = enrichGO(gene = GO,OrgDb = anas_org,keyType ="GO",ont = "MF",pvalueCutoff =0.05 ,qvalueCutoff = 0.05)这两步的意思是：
1. 首先通过orgdb进行更名，统一gene_symbol的名称。orgdb的本质就是一个巨大的命名表，里面是一个基因在各种数据库中一共24个各个类型的名称。
> columns(anas_org)[1] "ACCNUM" "ALIAS" "CHR" "ENSEMBL" "ENTREZID" "EVIDENCE" "EVIDENCEALL"[8] "GENENAME" "GID" "GO" "GOALL" "ONTOLOGY" "ONTOLOGYALL" "PMID"
[15] "REFSEQ" "SYMBOL" 2.org.db同时提供了索引背景基因的作用，因为GO富集作为超几何分布检验，无非是个背景基因和前景基因数量进行的一个计算。
问题就是出现在这里，这个非模式生物的org.db仅仅在我所导入的全物种基因中识别了几百个，仅仅是一个交集，这就很尴尬了。在我检查了自己之后，认定这个org.db数据库有点问题。（当然各位大佬要是发现我的步骤哪里不对的话请帮忙指出啊TAT）
2.通过uniprot自定义背景基因
当然，首先推荐一下刘小泽老师的博客，讲自己构建org.db的过程
https://www.jianshu.com/p/9c9e97167377/
刘小泽老师的步骤要下载并比对大量的基因来确定ID，但是我的数据量太大网络也太卡，一个简单的NR数据库一周都搞不定（就是下不完整怎么破）。
这里我们可以绕过org.db，方法是通过uniprot数据库的注释
1.首先，通过分类学找到自己需要的物种
2.点击Browse all浏览所有的基因信息
3.点击Customize columns来调整显示的信息
4.然后选择上GO注释
5.这下就可以看到所有基因的GO ID了
6.之后我们download下所有的信息（以表格的方式可以得到GO注释，别的好像不行）
7.打开R，开始编辑
接下来就是个命名的游戏了
整体思路：
通过enricher这个函数来解决这个问题
enricher函数需要输入三个数据：差异基因，名称转换表，描述信息
1.差异基因：一条基因id
2.名称转换表：两列数据框，第一列是基因id，第二列是对应的goid，只要保证大于差异基因列表就行
3.描述表：三列数据框：第一列是go id，第二列是基因描述，第三列是go分类（MF，CC，BP），是全部背景基因的集合。
注：需要实现根据go分类将名称转换表和描述表进行分类，因为enricher相当于一个做超几何分布的函数，他的运行逻辑就是：读入gene_id——转化goid——描述表中有多少基因富集到某个功能？——差异基因中有多少基因富集到某个功能？——超几何分布
根据这个逻辑，他是不会对goid背后的go分类进行解析的，所以需要事先通过go分类将uni2go和description来分成BP组，CC组，MF组
具体代码如下：
library(tidyr)
library(clusterProfiler)
library(stringr)
###########1.准备gene_symbol###########gene_symbol <- read.csv('C:/Users/gah20/Desktop/xxx/uniprot-txid8835_CDS_100979.csv',header = F)
gene_symbol<-gene_symbol$V17###########2.准备uni2go###########
uni2go<-read.csv('C:/Users/gah20/Desktop/xxx/my_data.csv',header = T)
View(head(uni2go,10))
colnames(uni2go)<-c('gene_id','ID')
# length(unique(uni2go$gene_id))
#>[1] 100979
uni2go <-uni2go %>% as_tibble() %>% separate_rows(ID, sep = ";")
remove_blank<-function(x){str_replace_all(x," ", "")}#这里是发现有些GO前面是带有空格，删除空格
uni2go$ID<-lapply(uni2go$ID,remove_blank)###########3.准备description文件########
dp<-read.csv('C:/Users/gah20/Desktop/xxx/description.csv',header = T)
dp<-unlist(strsplit(dp$锘縂ene.Ontology..GO.,split = ';'))#这个下下来，GOID那行的名称是乱码，就将乱就乱了，总之目的是把go_id1;go_id2;go_id3...这样的格式拆开
dp<-as.data.frame(dp)
dp$GO<-str_extract(dp$dp,"(?<=\\[).*?(?=\\])")
dp$description<-str_extract(dp$dp,".*?(?=\\[)")
dp<-dp[,-1]###########4.根据BP、MF、CC分割description与MF############
go_class<-read.csv('C:/Users/gah20/Desktop/xxx/description.csv',header = T)
View(head(go_class,10))
#uniprot数据库上下载的数据，是将MF，BP，CC分成了三列，然后把Go_id中的内容分配到其中，所以要用这种方法注释
MF<-go_class$Gene.Ontology..molecular.function.
MF<-MF[which(MF!="")]
MF <-MF %>% as.matrix() %>% as_tibble() %>%separate_rows(V1, sep = ";")BP<-go_class$Gene.Ontology..biological.process.
BP<-BP[which(BP!="")]
BP <-BP %>% as.matrix() %>% as_tibble() %>%separate_rows(V1, sep = ";")CC<-go_class$Gene.Ontology..cellular.component.
CC<-CC[which(CC!="")]
CC <-CC %>% as.matrix() %>% as_tibble() %>%separate_rows(V1, sep = ";")
MF<-str_extract(MF$V1,"(?<=\\[)GO.*?(?=])")%>% as.data.frame()
BP<-str_extract(BP$V1,"(?<=\\[)GO.*?(?=])")%>% as.data.frame()
CC<-str_extract(CC$V1,"(?<=\\[)GO.*?(?=])")%>% as.data.frame()colnames(MF)<-c("GO")
colnames(BP)<-c("GO")
colnames(CC)<-c("GO")MF2<-subset(dp,GO%in%MF$GO)
BP2<-subset(dp,GO%in%BP$GO)
CC2<-subset(dp,GO%in%CC$GO)uni2go_MF<-subset(uni2go,ID%in%MF2$GO)
uni2go_BP<-subset(uni2go,ID%in%BP2$GO)
uni2go_CC<-subset(uni2go,ID%in%CC2$GO)###########4.运行############
Go_rich_MF <- enricher(gene = gene_symbol,TERM2GENE = uni2go_MF[c('ID','gene_id')],TERM2NAME = MF2[c('GO','description')],pvalueCutoff = 0.05,pAdjustMethod = 'BH',qvalueCutoff = 0.2,maxGSSize = 200)Go_rich_BP <- enricher(gene = gene_symbol,TERM2GENE = uni2go_BP[c('ID','gene_id')],TERM2NAME = BP2[c('GO','description')],pvalueCutoff = 0.05,pAdjustMethod = 'BH',qvalueCutoff = 0.2,maxGSSize = 200)Go_rich_CC <- enricher(gene = gene_symbol,TERM2GENE = uni2go_CC[c('ID','gene_id')],TERM2NAME = CC2[c('GO','description')],pvalueCutoff = 0.05,pAdjustMethod = 'BH',qvalueCutoff = 0.2,maxGSSize = 200)相关的注释直接写进代码里面了
结语：本人也是处在一个摸索的阶段，各位路过的大佬有什么指点，请不吝赐教~
自定义背景基因非模式生物Go富集：使用uniprot数据自行建库，不需要orgdb相关推荐
- 基于对WalkGIS 的地形数据建库处理模式的探讨 来自建筑工程与设计 2015（15）
[摘要] 以1∶ 500 基础地理信息数据加工整理模式为研究对象,对数据整理的技术要求和整理方法进行了清晰的阐述,总结出一套基于WalkGIS 平台的基础地理信息数据库建库数据加工整理模式,为数据整理 ...
- Cell：康奈尔大学郭春君组开发针对非模式肠道细菌的基因编辑工具
北京时间2022年1月20日凌晨0时,美国康奈尔大学威尔康奈尔医学院郭春君(Chun-Jun Guo)研究组在<细胞>(Cell)上在线发表题为"Genetic manipula ...
- 模式生物，鸡--不止于大餐，还可用于基因编辑！- MedChemExpress
鸡 (学名:Gallus gallus domesticus) 是最常见的家禽之一,是全球肉类和蛋白质的重要来源.鸡体型小.易繁殖.孵化期短 (仅 21 天),性状优良的蛋鸡甚至在一年中的 316 天 ...
- MFC—对话框程序—模式对话框与非模式对话框
一.根据主窗口类型,MFC软件工程可以分为以下几种架构模型: 1.SDI(Single Document Interface):单文档界面,一个主框架窗口下只能编辑一份文档. 例如:记事本和画笔等. ...
- 文献集锦 | 非因生物空间多组学技术在头颈部肿瘤中的研究策略
头颈部鳞状细胞癌(HNSCC)是全球第七大癌症病因,是一种异质性恶性肿瘤,起源于上呼吸道,尤其是鳞状粘膜线.唇部.口腔和鼻腔.鼻窦.喉.鼻咽.口咽和下咽是HNSCC的受累部位.利用空间组学分析平台深入 ...
- 一个被遗忘的ccflow工作流引擎自定义表单开发模式
定义概述:一个已经做好的表单需要绑定到节点上 , 该文章在驰骋工作流引擎流程引擎设计器中. 自定义表单工作模式:流程控制按钮区域是ccflow来完成,表单区域是放在控制区域下面的框架里,如下图所示. ...
- Python库：wordcloud库介绍、政府工作报告词云、自定义背景词云
一.wordcloud库 二.使用wordcloud库 注:库名wordcloud全部是小写,而WordCloud对象W和C大写 简单说,绘制一个词云有三步: 第一.生成词云对象WordCloud,并 ...
- java 自定义形状按钮_制作自定义背景Button按钮、自定义形状Button的全攻略
在Android开发应用中,默认的Button是由系统渲染和管理大小的.而我们看到的成功的移动应用,都是有着酷炫的外观和使用体验的.因此,我们在开发产品的时候,需要对默认按钮进行美化.在本篇里,笔者结 ...
- pythonGUI实现照片或证件照迅速更换自定义背景底色
有一天,一位朋友发来威信吐槽道:唉,刚换了工作,又得换工作照,去拍照太麻烦,去PS也得花个个把小时还不一定修的好,去美图秀秀吧,竟然换一张照片的背景底色还要收我5.9元钱.于是老铁就问我有没有办法解决 ...
最新文章
- Nature：原来益生菌是这么搞定致病菌的
- java 中的内部类学习小记
- 基于MINA框架快速开发网络应用程序
- 大厂首发：kafka消费组订阅多个topic
- klib库下的kroundup32(二进制的四舍五入)算法
- java cocoon_Java-跳跃路线
- 创建Vue实例传入的options||Vue的生命周期
- 北京招聘 | 百度智能生活事业群组小度科技招聘对话系统算法实习生、工程师...
- 华为云 和 阿里云 跨服务器搭建Hadoop集群
- centos网络隔一段时间就断_“路由器隔一段时间就上不了网，断一下电又能用了，这是什么原因...
- 五、OpenStack安装Nova
- 物件捆绑 背包问题 动态规划 求解
- Drupal 7.31 SQL注入漏洞
- Tecplot 安装记录
- oracle节假日,oracle 产生节假日表
- .NET Core剪裁器背后的技术
- 批量剔除consul无效服务
- Cadence OrCAD Capture 在Excel中编辑所有元件属性然后导入设计图纸的方法
- python入门——Python的两种编程语言：交互式和文件式
- 基于阿里云IOT Studio和STM32的电机远程监测设计