REACTOME 是一个开源、放开、人工整理并经过同行评审的通路数据库。目标是为各个通路提供可视化、可解释的分析，以支持基础研究、临床研究、基因组分析、建模、系统生物学等。
和我们所熟知的GO/KEGG类似，Reactome这些年也越来越的的在文章中出现，截止至笔者写本文的时间，数据收录信息如下：
不同于KEGG这种收费的数据库，Reactome的使用非常的简单，现在让我们一步步的开始搭建分析工具吧！
首选需要安装以下包：clusterProfiler、ReactomePA、stringr、ggplot2，这些包均可以使用 BiocManager::install()或install.package()命令安装这里就不再赘述。下面直接贴出代码供大家使用。
suppressMessages(library(clusterProfiler))
suppressMessages(library(ReactomePA))
suppressMessages(library(stringr))
suppressMessages(library(ggplot2))
args <- commandArgs(T)
ref <- args[1]
genelist <- args[2]
outdir <- args[3]
top <- args[4]
GetGeneType<-function(df,filename){
if(str_detect(df[1,1],"EN")){
if( str_detect(df[1,1],"T") )
return("ENSEMBLTRANS")
else
return("ENSEMBL")
}else if(is.na(as.integer(df[1,1])) == F){
return("ENTREZID")
}else if(str_detect(df[1,1],"AT")){
return("TAIR")
}
else{
return("SYMBOL")
}
}
sample_name <- strsplit( basename(genelist) ,split=".",fixed=TRUE)[[1]][1]
df <- read.table(file=genelist,sep='\t',header=F)
setwd(outdir)
input_type = GetGeneType(df,genelist)
if(ref=='hs'){
suppressMessages(library(org.Hs.eg.db))
Db <- org.Hs.eg.db
organism <- 'hsa'
ref <- 'human'
}else if(ref=='mm'){
suppressMessages(library(org.Mm.eg.db))
Db <- org.Mm.eg.db
organism <- 'mmu'
ref <- 'mouse'
}else if( ref =='rn'){
suppressMessages(library(org.Rn.eg.db))
Db <- org.Rn.eg.db
organism <- 'rno'
ref <- 'rat'
}
bitrdf <- bitr(df[,1],fromType=input_type,toType="ENTREZID",OrgDb = Db)
Reactome <- as.data.frame(enrichPathway(gene=bitrdf[,'ENTREZID'] , pvalueCutoff = 1, readable=TRUE, organism = ref))
write.table(Reactome,paste0(sample_name,".Reactome.xls"),sep="\t",quote=F,row.names=F)
sigdf <- subset(Reactome, pvalue < 0.05)
write.table(sigdf,paste0(sample_name,".SigReactome.xls"),sep="\t",quote=F,row.names=F)
library(tidyverse)
plotdf <- separate(data = sigdf, col = BgRatio, into = c("BgCount", "BgAll"), sep = "/")
plotdf$BgCount = as.integer(plotdf$BgCount)
plotdf$RichFactor <- plotdf$Count / plotdf$BgCount
plotdf <- plotdf[order(plotdf$pvalue),]
plotdf <- head(plotdf,top)
p <- ggplot(data = plotdf,mapping = aes(x = RichFactor,y = reorder(Description,RichFactor)))+
geom_point(aes(color= pvalue,size = Count)) +
scale_colour_gradient(low = "red", high = "blue") +
theme_bw()+
labs(title = paste('Top',nrow(plotdf),'of Pathway Enrichment'),
x = 'Rich factor',
y = 'Pathway')+ scale_size("GeneNumber")+theme(plot.title = element_text(hjust = 0.5))
ggsave(paste0(sample_name,".bubble.png"), width = 8.5)
ggsave(paste0(sample_name,".bubble.pdf"), width = 8.5)
安装好对应的包后，将上述代码保存为Reactome.r使用命令：
Rscript Reactome.r hs gene.glist ./ 20 即可进行分析，其中
hs：物种信息，如人hs，小鼠mm，大鼠rn
gene.glist：基因列表一行一个基因，并去重，可以是ENSEMBL ID，ENTREZID, Symbol等，这里使用了bitrdf()将所有输入全部转为ENTREZID
./：输出路径，请先建立好该目录，./代表当前目录
20：绘制TopN的气泡图
下面是实际演示：
输入文件如下：
gene.glist文件
运行
Rscript Reactome.r hs gene.glist ./ 20
输出
输出文件
*.bubble.png/pdf是气泡图
气泡图
*.Reactome.xls和 *,SigReactome.xls是输出结果，前者是所有的富集结果后者是P<0.05的结果。
Reactome富集结果
这样我们就完成了整个Reactome富集分析了，心动了吗？那快加入你的文章当中吧！
