如何利用clusterProfiler获取最新的KEGG和基因对应关系

本文是洲更写的，文中说到两步，其实只需要这一步，然后就可以用enricher或GSEA函数进行分析，也就是通用的分析函数，不管是什么注释，不管是什么物种，都可以支持。

这样做的好处也有一些

• 分析时不用联网，可以搞个local版本

• 分析具有可重复性，实时爬网络数据，毕竟有时候服务器端有更新的话，就会导致有一些变化。

第二步所谓需要序列，其实是想本地blast，自己注释，这种一般是和KO数据库进行比对。等候洲更的更新，看他后面怎么玩。

---

Y叔的clusterProfiler的一大优点就是能够利用最新的KEGG数据库，而不是停留在最后一个公开版的KEGG数据库(2011-5-15).

大部分情况下，大家都是直接用enrichKEGG()或者gseKEGG()完成富集分析，但是我最近想到，我其实可以利用这个功能反向建立一个注释用的数据库，只需要两步

第一步: 获取KEGG编号对应的基因编号
第二步: 根据基因编号获取序列信息

这里只讲第一步，如何用Y叔clusterProfiler获取最新的KEGG和基因的对应关系， 我们以人类为例。
先用download_KEGG下载给定物种的KEGG数据库

hsa_kegg <- clusterProfiler::download_KEGG('hsa')


这一步得到的是一个列表，这个列表有两个成员，

names(hsa_kegg)
[1] 'KEGGPATHID2EXTID' 'KEGGPATHID2NAME' 

一个是KEGG的通路编号和基因编号的关系，另一个是KEGG通路编号和名字的关系

我们可以将其进行合并

PATH2ID <- hsa_kegg$KEGGPATHID2EXTID
PATH2NAME <- hsa_kegg$KEGGPATHID2NAME
PATH_ID_NAME <- merge(PATH2ID, PATH2NAME, by='from')
colnames(PATH_ID_NAME) <- c('KEGGID', 'ENTREZID', 'DESCRPTION')


最后的PATH_ID_NAME的表格信息如下

保存到本地

write.table(PATH_ID_NAME, 'HSA_KEGG.txt', sep='\t')

拓展：如何在此基础上增加ENSEMBL的编号？这可以使用biomaRt进行ID转换

library(biomaRt)

mart <- useDataset('hsapiens_gene_ensembl', useMart('ensembl'))
entrezgene <- PATH_ID_NAME$ENTREZID
# This step need some time
ensembl_gene_id<- getBM(attributes=c('ensembl_gene_id', 'entrezgene'),
                  filters = 'entrezgene',
                       values=entrezgene , mart= mart)


于是我们得到了ENSEMBL的基因编号，通过merge就可以将其添加到之前的数据框中

PATH_ID_NAME <- merge(PATH_ID_NAME, ensembl_gene_id, by.x= 'ENTREZID',by.y= 'entrezgene')

---