發(fā)現(xiàn)感興趣的KEGG ID居然不在KEGG.db包里面，比如：
hsa05034 Alcoholism
hsa05030 Cocaine addiction
導(dǎo)致下面的代碼失效：
library(KEGG.db)
ls("package:KEGG.db")
cellcycle_genes=KEGGPATHID2EXTID[['hsa04110']]
cytokine_genes=KEGGPATHID2EXTID[['hsa04060']]
KEGGPATHID2EXTID[['hsa05034']]
KEGGPATHID2EXTID[['hsa05030']]
其實，即使不失效，也不能使用這種方法來找屬于某個pathway的基因集合了，因為這個R包以及很多年沒有更新了。
搜索了一下，發(fā)現(xiàn)KEGG數(shù)據(jù)庫的rest API，比如
http://rest.kegg.jp/get/hsa05034 (點擊閱讀原文可以直達)
本來準備讀入到R里面，然后自己解析，發(fā)現(xiàn)其實已經(jīng)有了R包：
library(KEGGREST)
listDatabases()
## ----get_organisms------------------------------------------------------------
org <- keggList("organism")
head(org)
keggGet('hsa05034')
gs <- keggGet('hsa05034')
gs[[1]]$GENE
genes <- unlist(lapply(gs[[1]]$GENE,function(x) strsplit(x,';')[[1]][1]))
genes[1:length(genes)%%2 ==0]
當然了，這個R包的功能不止如此：https://bioconductor.org/packages/release/bioc/html/KEGGREST.html 我就不多演示了，感興趣的朋友去探索一下。
pathway gif動畫版 可視化 ：https://github.com/ajmazurie/kegg-animate-pathway
ensembl2symbol <- function(genes){
library(org.Hs.eg.db)
eg2symbol=toTable(org.Hs.egSYMBOL)
eg2ensembl=toTable(org.Hs.egENSEMBL)
#genes=sample(eg2ensembl$ensembl_id,100)
tmp=merge(eg2symbol,eg2ensembl,by='gene_id')
return(tmp[tmp$ensembl_id %in% genes,])
}
genes=sample(eg2ensembl$ensembl_id,100)
ensembl2symbol(genes)
genes=sample(eg2ensembl$ensembl_id,100)
ensembl2symbol(genes)
genes=sample(eg2ensembl$ensembl_id,100)
ensembl2symbol(genes)
就當是我當年參加生信技能樹直播編程活動提交的作業(yè)吧！
歷史題目：
生物信息學(xué)技能面試題(第1題)-人類基因組的外顯子區(qū)域到底有多長
生物信息學(xué)技能面試題(第2題)-探索人類基因組序列
生物信息學(xué)技能面試題(第3題)-探索人類基因組注釋文件
生物信息學(xué)技能面試題(第4題)-多個同樣的行列式文件合并起來
生物信息學(xué)技能面試題(第5題)-根據(jù)GTF畫基因的多個轉(zhuǎn)錄本結(jié)構(gòu)
生物信息學(xué)技能面試題(第6題)-下載最新版的KEGG信息，并且解析好
生信編程直播第9題-根據(jù)指定染色體及坐標得到參考堿基
生信編程直播第10題：根據(jù)指定染色體及坐標得到位置信息
聯(lián)系客服
