国产一级a片免费看高清,亚洲熟女中文字幕在线视频,黄三级高清在线播放,免费黄色视频在线看

打開APP
userphoto
未登錄

開通VIP,暢享免費電子書等14項超值服

開通VIP
KEGG數(shù)據(jù)庫倒閉了嗎

最近,有粉絲運行了我以前的數(shù)據(jù)挖掘成套代碼里面的 run_kegg 函數(shù),如下所示:

library(clusterProfiler) 
run_kegg(gene_up,gene_down,pro='comp1')

出現(xiàn)了如下所示的報錯:

Reading KEGG annotation online:

fail to download KEGG data...
Error in download.KEGG.Path(species) : 
  'species' should be one of organisms listed in 'http://www.genome.jp/kegg/catalog/org_list.html'...
In addition: Warning message:
In utils::download.file(url, quiet = TRUE, method = method, ...) :
  URL 'https://rest.kegg.jp/link/mmu/pathway': status was 'Failure when receiving data from the peer'
Called from: download.KEGG.Path(species)

然后就找我,以為是我們的標準代碼有問題,實際上我的 run_kegg 函數(shù)僅僅是包裝了 Y叔的 clusterProfiler包而已 ,實際上里面沒有啥玄機,如下所示:

## KEGG pathway analysis
### 做KEGG數(shù)據(jù)集超幾何分布檢驗分析,重點在結(jié)果的可視化及生物學意義的理解。
run_kegg <- function(gene_up,gene_down,pro='test'){
  library(ggplot2)
  gene_up=unique(gene_up)
  gene_down=unique(gene_down)
  gene_diff=unique(c(gene_up,gene_down))
  ###   over-representation test
  # 下面把3個基因集分開做超幾何分布檢驗
  # 首先是上調(diào)基因集。
  kk.up <- enrichKEGG(gene   = gene_up,
 organism  = 'mmu',
 #universe  = gene_all,
 pvalueCutoff = 0.9,
 qvalueCutoff =0.9)
  head(kk.up)[,1:6]
  kk=kk.up
  dotplot(kk)
  kk=DOSE::setReadable(kk, OrgDb='org.Mm.eg.db',keyType='ENTREZID')
  write.csv(kk@result,paste0(pro,'_kk.up.csv'))
  
  # 首先是下調(diào)基因集。
  kk.down <- enrichKEGG(gene   =  gene_down,
   organism  = 'mmu',
   #universe  = gene_all,
   pvalueCutoff = 0.9,
   qvalueCutoff =0.9)
  head(kk.down)[,1:6]
  kk=kk.down
  dotplot(kk)
  kk=DOSE::setReadable(kk, OrgDb='org.Mm.eg.db',keyType='ENTREZID')
  write.csv(kk@result,paste0(pro,'_kk.down.csv'))
  
  # 最后是上下調(diào)合并后的基因集。
  kk.diff <- enrichKEGG(gene   = gene_diff,
   organism  = 'mmu',
   pvalueCutoff = 0.05)
  head(kk.diff)[,1:6]
  kk=kk.diff
  dotplot(kk)
  kk=DOSE::setReadable(kk, OrgDb='org.Mm.eg.db',keyType='ENTREZID')
  write.csv(kk@result,paste0(pro,'_kk.diff.csv')) 
}

就是針對上下調(diào)基因,以及其合并的差異基因,分別獨立走 enrichKEGG 函數(shù)而已,它來自于 Y叔的 clusterProfiler包。

實際上,我們根據(jù)前面的報錯,進入  https://www.genome.jp/kegg/catalog/org_list.html 可以很清楚的看到,我們的人和鼠,都是在里面的:

hsa   Homo sapiens (human)
mmu   Mus musculus (house mouse)

所以不可能是物種問題,而且這個kegg數(shù)據(jù)庫的官網(wǎng)也可以訪問,那么合理的推測,是不是Y叔的 clusterProfiler包有問題呢?其實仔細看報錯提示是 download.KEGG.Path 函數(shù)問題,也是Y叔的包里面的。這個download.KEGG.Path 函數(shù)的代碼是:

 keggpathid2extid.df <- kegg_link(species, "pathway")
 if (is.null(keggpathid2extid.df)) 
  stop("'species' should be one of organisms listed in 'http://www.genome.jp/kegg/catalog/org_list.html'...")
 keggpathid2extid.df[, 1] %<>% gsub("[^:]+:"""
  .)
 keggpathid2extid.df[, 2] %<>% gsub("[^:]+:"""
  .)
 keggpathid2name.df <- kegg_list("pathway")
 keggpathid2name.df[, 1] %<>% gsub("path:map", species, 
  .)
 keggpathid2extid.df <- keggpathid2extid.df[keggpathid2extid.df[, 
  1] %in% keggpathid2name.df[, 1], ]
 return(list(KEGGPATHID2EXTID = keggpathid2extid.df, KEGGPATHID2NAME = keggpathid2name.df))

也就是說在kegg_link函數(shù)就開始報錯了,繼續(xù)看kegg_link的代碼:

url <- paste0("http://rest.kegg.jp/link/", target_db, 
  "/", source_db, collapse = "")
 kegg_rest(url)

實際上, 它就根據(jù)我們的物種,拼湊了一個 網(wǎng)頁,比如:https://rest.kegg.jp/link/hsa/pathway,可以使用瀏覽器打開,能看到里面的通路的id和基因id的對應(yīng)關(guān)系。所以真正的問題是 kegg_rest 這個函數(shù),沒有辦法像瀏覽器那樣直接獲取網(wǎng)頁里面的內(nèi)容。我們繼續(xù)看kegg_rest 這個函數(shù)的源代碼:

 dl <- mydownload(rest_url, destfile = f)
 if (is.null(dl)) {
  message("fail to download KEGG data...")
  return(NULL)
 }

真氣人,就是一個套娃, 里面又是mydownload函數(shù),那我們繼續(xù)看mydownload函數(shù)的源代碼 :

if (is.null(method)) 
  method <- getOption("clusterProfiler.download.method")
 if (!is.null(method) && method != "auto") {
  dl <- tryCatch(utils::download.file(url, quiet = TRUE
   method = method, ...), error = function(e) NULL)
 }
 else {
  dl <- tryCatch(downloader::download(url, quiet = TRUE
   ...), error = function(e) NULL)
 }
 return(dl)

有意思的是,我看了看,它里面的兩個下載函數(shù)(utils::download.file和  downloader::download)我都測試了,是可以獨立下載的:

utils::download.file( "http://rest.kegg.jp/link/hsa/pathway" ,
tempfile()  )
#  downloader::download 也可以

trying URL 'http://rest.kegg.jp/link/hsa/pathway'
Content type 'text/plain; charset=utf-8' length unknown
downloaded 807 KB

但是這兩個下載函數(shù)(utils::download.file和  downloader::download),我是使用默認參數(shù)進行下載,都沒有設(shè)置協(xié)議,我看了看函數(shù)里面的協(xié)議是:

 getOption("clusterProfiler.download.method")
[1] "libcurl"

也就是說,它們兩個下載函數(shù)(utils::download.file和  downloader::download)使用默認的方法,也就是 'auto'  是可以去訪問我們的瀏覽器可以訪問的:https://rest.kegg.jp/link/hsa/pathway,的內(nèi)容。但是如果給它們 "libcurl"的協(xié)議,就會失敗。

所以,最簡單的解決方案就是修改這個默認的協(xié)議即可,我給一個解決方案,就是:

install.packages('R.utils')
R.utils::setOption( "clusterProfiler.download.method",'auto' )

果然,接下來的富集分析就非常絲滑了。

出圖如下所示:

 

KEGG數(shù)據(jù)庫沒有倒閉, Y叔的 clusterProfiler包也問題不大,我的一個 run_kegg 函數(shù)更不可能有問題。僅僅是因為R語言里面的下載文件的函數(shù)的協(xié)議需要注意,這兩個函數(shù)兩個下載函數(shù)(utils::download.file和  downloader::download),都太底層了。初學者確實不容易找到問題所在。

寫在文末

我在《生信技能樹》,《生信菜鳥團》,《單細胞天地》的大量推文教程里面共享的代碼都是復制粘貼即可使用的, 有任何疑問歡迎留言討論,也可以發(fā)郵件給我,詳細描述你遇到的困難的前因后果給我,我的郵箱地址是 jmzeng1314@163.com

如果你確實覺得我的教程對你的科研課題有幫助,讓你茅塞頓開,或者說你的課題大量使用我的技能,煩請日后在發(fā)表自己的成果的時候,加上一個簡短的致謝,如下所示:

We thank Dr.Jianming Zeng(University of Macau), and all the members of his bioinformatics team, biotrainee, for generously sharing their experience and codes.

十年后我環(huán)游世界各地的高校以及科研院所(當然包括中國大陸)的時候,如果有這樣的情誼,我會優(yōu)先見你。

本站僅提供存儲服務(wù),所有內(nèi)容均由用戶發(fā)布,如發(fā)現(xiàn)有害或侵權(quán)內(nèi)容,請點擊舉報
打開APP,閱讀全文并永久保存 查看更多類似文章
猜你喜歡
類似文章
clusterProfiler——GO和KEGG分析之R代碼
KEGG數(shù)據(jù)本地化,再也不用擔心網(wǎng)絡(luò)問題了
GEO數(shù)據(jù)挖掘流程——代碼版(方便抄襲)
KEGG數(shù)據(jù)庫的rest API(附帶R語言小技巧)
復現(xiàn)Nature圖表:GSEA分析及可視化包裝函數(shù)
得到差異分析之后進行功能富集分析
更多類似文章 >>
生活服務(wù)
分享 收藏 導長圖 關(guān)注 下載文章
綁定賬號成功
后續(xù)可登錄賬號暢享VIP特權(quán)!
如果VIP功能使用有故障,
可點擊這里聯(lián)系客服!

聯(lián)系客服