方法一
在MSigDB下载,网址GSEA | MSigDB | MSigDB Collections。下载到 gmt
格式文件使用 GSEABase
包的 getGmt
函数读取。这个下的可能有点旧,目前(2019.11.27)版本是v7.0共186条KEGG通路,显然比目前KEGG网站收录的少。
方法二
从注释包 org.Hs.eg.db
里提取。目前 org.Hs.eg.db
版本是 3.8.2 从里面提取到229条通路。从注释包提取信息方法参见我上一篇文章《AnnotationDbi 使用(以 org.Hs.eg.db 为例)》。
方法三
用R包 KEGGREST
批量检索下载。 KEGGREST
是个通过KEGG API在线获取KEGG数据库信息的工具,因此相对来说能够获得更新更全的数据。通过3步就能获取到所有的KEGG通路,1获取KEGG通路列表,2获取每个通路基因,3整理数据。大致代码如下,详细解释看代码注释。
library(KEGGREST, quietly = TRUE)
library(tidyverse, quietly = TRUE)
# 返回信息很长,只取基因symbol.根据自己需要调整
symbolOnly <- function(x){
items <- strsplit(x, ";", fixed = TRUE) %>% unlist()
return(items[1])
}
# keggGet(x)[[1]]$GENE 数据基因名是个向量,其中奇数位置是 entrezgene_id 偶数位置是 symbol
geneEntrez <- function(x){
geneList <- keggGet(x)[[1]]$GENE
if(!is.null(geneList)){
listLength <- length(geneList)
entrezList <- geneList[seq.int(from = 1, by = 2, length.out = listLength/2)]
entrez <- stringr::str_c(entrezList, collapse = ",")
return(entrez)
}else{
return(NA)
}
}
# keggGet(x)[[1]]$GENE 数据基因名是个向量,其中奇数位置是 entrezgene_id 偶数位置是 symbol
geneSymbol <- function(x){
geneList <- keggGet(x)[[1]]$GENE
if(!is.null(geneList)){
listLength <- length(geneList)
symbolList <- geneList[seq.int(from = 2, by = 2, length.out = listLength/2)] %>% map_chr(symbolOnly)
symbol <- stringr::str_c(symbolList, collapse = ",")
return(symbol)
}else{
return("")
}
}
# 取得 hsaxxxxx 这种通路ID
pathwayID <- function(x){
items <- strsplit(x, ":", fixed = TRUE) %>% unlist()
return(items[2])
}
# 建议从这里开始读脚本。建议自己在交互模式下试一下用到的KEGGREST函数,看看返回数据的结构。
# 这是第一步,取得所有的KEGG通路列表
hsaList <- keggList("pathway", "hsa")
IDList <- names(hsaList) %>% map_chr(pathwayID)
# 将通路ID和通路名放在一个表格(tibble)里
hsaPathway <- tibble::tibble(pathway_id=IDList, pathway_name=hsaList)
head(hsaPathway, n=3) %>% print()
# 用前面定义函数,获得每个通路的基因,然后也放在表格里
pathwayFull <- hsaPathway %>% dplyr::mutate(entrezgene_id=map_chr(pathway_id, geneEntrez), hgnc_symbol=map_chr(pathway_id, geneSymbol))
# 保存数据
write_tsv(pathwayFull, path="KEGGREST.tsv")
dim(pathwayFull) %>% print()
# 会有通路没有基因,我的话只需要有基因的,所以把无基因的移除
pathwayWithGene <- dplyr::filter(pathwayFull, !is.na(entrezgene_id) & hgnc_symbol != "")
write_tsv(pathwayWithGene, path="KEGGREST_WithGene.tsv")
dim(pathwayWithGene) %>% print()
网友评论