（未测试）数据挖掘 | KEGG通路与基因的对应关系

基于 KEGGREST 包：

p_load(KEGGREST)

以 hsa05215 (Prostate cancer) 通路为例。一切都要从 keggGet() 获取的全部信息开始：

infos = keggGet("hsa05215")

数据很全面，包括其相关联的通路、发表的相关文献、相关疾病、治疗药物、内部基因等等。

提取基因信息：

gene_info = matrix(infos[[1]]$GENE, ncol = 2, byrow = T)

如图所示，第一列为Gene ID，第二列为基因描述信息，以”;”分隔，前面即为Symbol，使用 stringr 包中的 str_extract() 函数可以很方便的提取出来：

p_load(stringr)str_extract(gene_info[,2], "^[^;]*")

共计97个基因。

拓展来看，如果想获得所有pathway对应的所有基因也是可行的，实现上并非针对每个hsa*编号重复执行上述过程，而是通过 keggLink() 函数：

kl = keggLink("pathway", "hsa")
path_gid = data.frame(pathway = str_extract(kl, "[^:]*$"), id = str_extract(names(kl), "[^:]*$")

如图所示，第一列为通路编号，第二列为Gene ID。结合 org.Hs.eg.db 包，我们还可以做更多，比如获取Symbol和Ensembl ID:

p_load(org.Hs.eg.db, tidyverse)
path_gene = path_gid %>%
    mutate_all(as.character) %>%
    mutate(symbol = mapIds(org.Hs.eg.db, id, "SYMBOL", "ENTREZID"),        
           ensembl = mapIds(org.Hs.eg.db, id, "ENSEMBL", "ENTREZID")
    )

再获取通路的描述信息：

paths = keggList("pathway", "hsa")
path_descriptions = data.frame(pathway = str_extract(names(paths), "[^:]*$"), description = str_match(paths, "(.*) - ")[,2])

合并整理下：

path_final = full_join(path_descriptions, path_gene, by = "pathway")

这个时候，通路和基因的对应关系就很明了了吧！得到这个表之后，我们至少还有三件事情可以做：

1、直接从中提取开题的hsa05215通路中的基因：

path_final %>% filter(pathway == "hsa05215") %>% dplyr::select(symbol)

也是97个基因，与之前的结果一致！

2、获取特定基因对应的通路，即哪些通路中包含着你要做的目标基因！

3、保存到文件中，下次可以直接离线使用！