1.安装

source("http://bioconductor.org/biocLite.R")
biocLite('biomaRt')
biocLite('clusterProfiler')

R 版本为3.5以上的需要以下命令

BiocManager::install('biomaRt')
BiocManager::install('clusterProfiler')

如果提示缺哪个包，则用以下命令安装：

install.packages('')#''内为缺失的包名

2.生物学ID转换

——做KEGG和GO分析之前需要对非NCBI ID的转换为NCBI ID，如山羊的Hoxc13基因的NCBI ID为102178245，人的Hoxc13基因的NCBI ID为3229。如果你在NCBI https://www.ncbi.nlm.nih.gov/的gene这一栏中输入102178245，那么它反馈回来的为山羊的Hoxc13基因，如果你输入3229，那么它反馈回来的为人的Hoxc13基因。

2.1 选择数据库

library(biomaRt)
mart <- useMart("ensembl","hsapiens_gene_ensembl")

——useMart第二个双引号内为需要使用的数据库，根据物种的不同，选择不同的数据库。如，人的选择hsapiens_gene_ensembl，鼠的选择mmusculus_gene_ensembl

2.2 选择输入ID类型

——然后，你需要知道输入的是什么类型的ID，即转换前的ID类型，如：

ENSG00000000003或ENMUSG000000003，属于类型为ensembl_gene_id；
ENST00000000233或ENMUST00000000233，属于类型为ensembl_transcript_id；
102178245，属于类型为entrezgene；
Hoxc13，属于类型为external_gene_name；
NM_000014，属于类型为refseq_mrna；
hsa-let-7a-1，属于类型为mirbase_id；

要想知道biomaRt支持哪些ID类型的输入，可以通过以下命令查看，共380种ID。

listFilters(mart)#输入的ID类型

2.3 选择输出ID类型

——要想知道biomaRt支持哪些ID类型的输出，可以通过以下命令查看，共支持2722种ID输出。

listAttributes(mart)

常用的：

Gene description：
Gene name：
Chromosome/scaffold:
Gene start (bp):
Gene end (bp):
KEGG Pathway and Enzyme:
NCBI gene:
PDB:
entrezgene_id,即NCBI gene ID
external_gene_name

2.4 导入文件和选定向量

gene <- read.csv(“data.csv”,sep=',',header=T)[,4]

因为需要转换的ID在数据集data.csv的第四列。

题外话

“data.csv”只是我自己的一个例子，如果你的为xlxs格式，可以另存为.csv格式，然后在通过read.csv函数导入；如果你的为txt格式，通过read.table函数导入，如：

gene <- read.table(“data.txt”,sep='\t',header=T)

其中双引号内的“data.txt”，有两种书写方法：
1，输入文件的完整路径（注意反斜线），如：“G:/Download/data.csv”
2，通过R菜单将工作目录切换到data.csv所在目录，“”内只写文件名即可。

gene # 查看结果为向量
#[1] NFKB1   CDKN1A  MMP16   KIT     Card10  Scube2  TRAF6   IRAK1   TLR4    UHRF1   PDGFRA  S100A12 ZNF117  MYO6   
#[15] CCDC6   IL1RAP  IL1RL2  HNRNPD  ZNRF3   IL6     PAX8    NOVA1   RARB    EGFR    ERBB4   SLC5A5  MALAT1  BSG    
#[29] YWHAZ   WEE1    TRPS1   MBNL2   KLF13   CDKN1A  LATS2   ERBB4   NR4A2   VEGFA   TGFBR2  RHOC    NFIB    CDK2   
#[43] CCNA1   TNFAIP1 DKK1    BTG1    LEFTY1  TXNIP   ATAD2  
#47 Levels: ATAD2 BSG BTG1 Card10 CCDC6 CCNA1 CDK2 CDKN1A DKK1 EGFR ERBB4 HNRNPD IL1RAP IL1RL2 IL6 IRAK1 KIT ... ZNRF3

2.5 开始转换

geneID <- getBM(attributes=c("entrezgene_id"),filters = "external_gene_name",values = gene, mart = mart)[,1]
geneID # 查看
#[1]      7126      8808      8900      6528      7048      2066 109504726     51351     26524     84133      7189
#[12]      1026     29775     29128     22943     57758      3569      4325      7227    378938      7422      5156
#[23]     29028     51621      4790     10628       389       694      3815      7849      7465      5915      3556
#[34]      8030      3184      4646      1956      4857      4781      7099      3654       682     10637      4929
#[45]      6283     10150      7534      1017

filters，指定输入ID的类型，
attributes指定输出ID的类型，
[,1]，取第一列，是为了把输出的ID从矩阵变为向量，因为clusterProfiler包的enrichKEGG函数需要输入为向量。
这次结果中，输入49个"external_gene_name"，输出48个"entrezgene"。

2.6 clusterProfiler使用

2.6.1 KEGG分析

library('clusterProfiler') 
kegg <- enrichKEGG(geneID, organism = "hsa", keyType = "kegg", use_internal_data = FALSE)

如果是使用clusterProfiler包的bitr函数进行ID转换，则可以这样做：

library('org.Hs.eg.db')
library('org.Mm.eg.db')#物种为鼠
ENTREZID<- bitr(gene, fromType = "SYMBOL", toType="ENTREZID",OrgDb = org.Hs.eg.db)
#Warning message:
#In bitr(gene, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = org.Hs.eg.db) :
#4.35% of input gene IDs are fail to map...

此外还有"ENSEMBL"。

geneID<-ENTREZID[,2]#取第二列
geneID
#[1] "4790"   "1026"   "4325"   "3815"   "7189"   "3654"   "7099"   "29128"  "5156"   "6283"   "51351"  "4646"   "8030"  
#[14] "3556"   "8808"   "3184"   "84133"  "3569"   "7849"   "4857"   "5915"   "1956"   "2066"   "6528"   "378938" "682"   
#[27] "7534"   "7465"   "7227"   "10150"  "51621"  "26524"  "4929"   "7422"   "7048"   "389"    "4781"   "1017"   "8900"  
#[40] "7126"   "22943"  "694"    "10637"  "10628"  "29028"
kegg <- enrichKEGG(geneID, organism = "hsa", keyType = "kegg", use_internal_data = TRUE)

可以看到clusterProfiler包输入49个 "SYMBOL",输出45个"ENTREZID"。所以有4.35% of input gene IDs are fail to map。"SYMBOL"和biomaRt包的"external_gene_name"，"ENTREZID"和biomaRt包的"entrezgene"其实是同一种类型，只不过两个包名字不同而已。
https://www.genome.jp/kegg/catalog/org_list.html查看支持的物种（organisms）。
人：hsa；
鼠：mmu。

2.5.2 GO分析

GO_BP <- enrichGO(geneID, OrgDb = org.Hs.eg.db, ont='BP',pAdjustMethod = 'BH', qvalueCutoff = 0.2, pvalueCutoff = 0.05,keyType = 'ENTREZID')
GO_CC <- enrichGO(geneID, OrgDb = org.Hs.eg.db, ont='CC',pAdjustMethod = 'BH', qvalueCutoff = 0.2, pvalueCutoff = 0.05,keyType = 'ENTREZID')
GO_MF <- enrichGO(geneID, OrgDb = org.Hs.eg.db, ont='MF',pAdjustMethod = 'BH', qvalueCutoff = 0.2, pvalueCutoff = 0.05,keyType = 'ENTREZID')

2.6 作图

dotplot(kegg,showCategory=20) #气泡图
barplot(kegg,showCategory=20,drop=T) #showCategory,显示条目的数目
barplot(GO_BP,showCategory=20,drop=T)
barplot(GO_CC,showCategory=20,drop=T)
barplot(GO_MF,showCategory=20,drop=T)

从biomaRt到clusterProfiler