你只想做ID转换却不知道为什么要转换

Posted on by

最近咱们《生信技能树》的VIP答疑群，有这样的提问：

image-20210520195748492

我觉得很有代表性，很多人仅仅是学了个皮毛，知道是需要进行ID转换，也能够运行代码。但是却搞不懂，不理解自己为什么进行ID转换，以及做了ID转换后是为了方便什么样的后续分析！

ID转换后我们才认识基因功能

在教程《研究最热门的基因是什么》里面，我提到了普通的词云展现的基因是 Entrez Gene，绝大部分人是没办法一目了然的去认识它，比如这个时候虽然知道研究最多的是7157和 7124 这两个基因，但是一般人是没办法认识这两个数字背后的基因啦。

所以，我们需要理所当然的做一个ID转换，如下所示的代码：

library(org.Hs.eg.db)
ids=toTable(org.Hs.egSYMBOL)
head(ids)
tbs=merge(ids,tb,by='gene_id')
wordcloud(words = tbs$symbol, freq = tbs$Freq, min.freq = 1,
 max.words=200, random.order=FALSE, rot.per=0.35,
 colors=brewer.pal(8, "Dark2"))

重新出图很清晰的看到了最热门的基因就是TP53

ID转换才能与其它数据库进行对接

这就是我们常见的GO/KEGG的数据库注释分析，我们这个时候并不会关心具体基因的名字，因为哪怕是正常的基因symbol，如果是成百上千个基因给到你，人力也是不可能是去全部记忆。但是把他们注释到数据库，就方便看他们整体的功能啦！

在：差异分析要的是表达量矩阵，基因名字并不重要啊，我提到过其实这样的注释并不一定要进行ID转换，因为只需要ID是一以贯之的即可：

image-20210522181956933

通常情况下，我们是拿到了上下调基因集合，去做GO和KEGG数据库注释，如果是人类数据，下面是一个示例代码：

gene_up=rownames(deg[deg$avg_logFC>0,])
gene_down=rownames(deg[deg$avg_logFC < 0,])

library(org.Hs.eg.db)
#把SYMBOL改为ENTREZID，这里会损失一部分无法匹配到的
gene_up=as.character(na.omit(AnnotationDbi::select(org.Hs.eg.db,
 keys = gene_up,
 columns = 'ENTREZID',
 keytype = 'SYMBOL')[,2]))
gene_down=as.character(na.omit(AnnotationDbi::select(org.Hs.eg.db,
 keys = gene_down,
 columns = 'ENTREZID',
 keytype = 'SYMBOL')[,2]))

# 人家的包就是 entrez ID 设计
library(clusterProfiler) 
#对正相关基因进行富集分析
go <- enrichGO(gene_up, OrgDb = "org.Hs.eg.db", ont="all") 
library(ggplot2)
library(stringr)
barplot(go, split="ONTOLOGY")+ facet_grid(ONTOLOGY~., scale="free")+ 
 ggsave('gene_up_GO_all_barplot.png') 
go=DOSE::setReadable(go, OrgDb='org.Hs.eg.db',keyType='ENTREZID')
write.csv( go@result,file = 'gene_up_GO_all_barplot.csv')

#对负相关基因进行富集分析
go <- enrichGO(gene_down, OrgDb = "org.Hs.eg.db", ont="all") 
barplot(go, split="ONTOLOGY")+ facet_grid(ONTOLOGY~., scale="free")+ 
 ggsave('gene_down_GO_all_barplot.png') 
go=DOSE::setReadable(go, OrgDb='org.Hs.eg.db',keyType='ENTREZID')
write.csv( go@result,file = 'gene_down_GO_all_barplot.csv')

如果是其它物种，尤其是小众的，可能是得自己慢慢摸索了。

那么lncRNA为什么可以不进行ID转换呢？

如果你是需要做GO/KEGG的数据库注释分析，其实绝大部分lncRNA本来就是并不会出现在主流的GO/KEGG的数据库，你把它的ID变来变去都是一样的结果哈！

大量lncRNA的功能是未知的，但是它们主要是cis-regulators，所以可以根据它们临近的蛋白编码基因功能来近似推断，然后表达量的相关性也可以类推到。

  • 如果是根据位置关系推断，使用bedtools等工具！

  • 如果是感觉表达量的相关性，可以看数据库。比如杂志Cancer Medicine, 2020的文章《 Genome-wide DNA methylation analysis by MethylRad and the transcriptome profiles reveal the potential cancer-related lncRNAs in colon cancer》，在进行结直肠癌相关lncRNA的功能富集分析，就是采用LncRN2Target v2.0和StarBase分析与15个lncRNA共表达的蛋白编码基因，其中lncRNA HULC和ZNF667-AS1分别鉴定到28个、9个共表达的蛋白编码基因！

​

Comments are closed.