功能富集分析流程(一)

功能富集分析流程是怎樣的?

接下來,我們以KEGG數據庫為例,帶大家梳理下功能富集分析流程

1、教你怎么將KEGG數據庫本地化。

2、如何自己定義功能富集的背景基因集(以物種人類為例,選取一個平臺檢測的所有基因和全部的人類已知基因做背景基因有什么區別么?差異在哪里?)

3、如何自己寫代碼實現功能富集核心原理。

小編有話說

首先做一個基因集合的功能富集,必須搞清楚一件事情,基因功能富集和基因功能注釋的差別:前者是一個 gene set 是否共同參與完成了一個生物學過程,而功能注釋只是針對single gene來說的,這個基因可能參與的生物學過程,有一些什么功能,只需要看這個基因是否在某個通路中出現,一個簡單的map過程而已,沒有用到統計學模型,不需要顯著性計算。

而評價一個gene set是否共同參與完成了某個生物學過程,這有很多統計學模型和算法來評估這個事情,比如最常見的超幾何模型和GSEA功能富集的排秩原理。這里小編會以超幾何原理為例,告訴你怎么自己寫代碼實現這個功能富集過程,加深對功能富集的理解。

第一步 KEGG數據庫本地化

》獲取KEGG通路中所有通路列表

1.1?KEGG官網https://www.genome.jp/kegg/,首先我們看一下KEGG數據庫中的Pathway。

KEGG將pathway分為了七大類,每一大類下面有小類,大家可以自己去詳細看看是什么樣的。

1.?Metabolism?代謝 代謝
2.?Genetic?Information?Processing 遺傳信息處理
3.?Environmental?Information?Processing 環境信息處理
4.?Cellular?Processes 細胞過程
5.?Organismal?Systems 組織系統
6.?Human?Diseases 疾病
7.?Drug?Development 藥物

1.2?進入如下這個網址:

https://www.kegg.jp/kegg/rest/keggapi.html#list,點擊如圖所示鏈接得到KEGG所有人類通路的列表,截至20190418,收錄了330條通路。

通路列表如下,第一列是KEGG通路ID號,第二列是通路名稱。

path:hsa00010 Glycolysis / Gluconeogenesis – Homo sapiens (human)
path:hsa00020 Citrate cycle (TCA cycle) – Homo sapiens (human)
path:hsa00030 Pentose phosphate pathway – Homo sapiens (human)
path:hsa00040 Pentose and glucuronate interconversions – Homo sapiens (human)
path:hsa00051 Fructose and mannose metabolism – Homo sapiens (human)
path:hsa00052 Galactose metabolism – Homo sapiens (human)
path:hsa00053 Ascorbate and aldarate metabolism – Homo sapiens (human)
path:hsa00061 Fatty acid biosynthesis – Homo sapiens (human)
path:hsa00062 Fatty acid elongation – Homo sapiens (human)

1.3?獲取任意一個pathway下載網址

進入網址https://www.kegg.jp/kegg/kegg3a.html,點擊圖片中框框位置,得到通路下載鏈接:

https://www.kegg.jp/kegg-bin/download?entry=hsa00010&format=kgml

1.4? 根據上面下載連接的特點和1.2通路列表id生成所有的KEGG通路下載地址

https://www.kegg.jp/kegg-bin/download?entry=hsa00010&format=kgml

即替換這個網址中entry=hsa00010即可,R代碼如下:

#?生成某一通路的鏈接,最后用迅雷批量下載xml文件
rm(list=ls())
setwd(“/path”)

#?導入kegg文件(通路ID,通路名字)
kegg?<-?read.table(“KEGG_pathwayID.txt”,header=F,sep=“\t”,quote=“\””,stringsAsFactors=F)
dim(kegg)?

#?生成某一通路的鏈接(手動網頁復制)
original_URL?<-?“http://www.kegg.jp/kegg-bin/download?entry=hsa00010&format=kgml”??
replaced_pathwayId?<-?gsub(“path:”,“”,kegg[,1])

replaced_URL?<-?sapply(replaced_pathwayId,function(x){gsub(“hsa00010”,x,original_URL)})

hsa_keggURLs?<-?cbind(pathwayId=replaced_pathwayId,pathwayName=kegg[,2],pathwayURL=replaced_URL)

write.table(hsa_keggURLs,file=“hsa_keggURLs.txt”,sep=“\t”,row.names=F,quote=F)

然后使用wget 命令下載hsa_keggURLs.txt文件中所有通路的xml文件,或者使用迅雷下載也是可以噠!

本期功能富集分析流程先講到這里,下期我們講解xml文件的特點和如何根據通路的xml文件提取每個通路中注釋的geneid,并獲得對應關系。

最近文章
女皇之心试玩