基因功能注释以及富集分析

生信喵

一、注释相关包介绍
       clusterProfiler 是南方医科大学余光创（微信公众号 biobabble）写的一个 R包，可以用来做各种富集分析，如 GO、KEGG、DO（Disease Ontology
analysis）、Reactome pathway analysis 以及 GSEA 富集分析等。而除了富集分析，还可以非常方便的对富集分析结果进行可视化。只需要输入一个基因 ID的列表，就可以完成的分析，非常方便。

网址：

1. http://www.bioconductor.org/packages/release/bioc/html/clusterProfiler.html

复制代码

使用说明：

1. http://www.bioconductor.org/packages/release/bioc/vignettes/clusterProfiler/inst/doc/clusterProfiler.html

复制代码

1. http://yulab-smu.top/clusterProfiler-book/

复制代码

安装和使用：

1. library(dplyr)
   
2. library(tidyr)
   
3. library(DOSE)
   
4. library(GO.db)
   
5. library(org.Hs.eg.db)
   
6. library(GSEABase)
   
7. library(clusterProfiler)

复制代码

二、基因 ID 转换
       对于具有注释信息的物种，基因注释的过程就是各种 ID 转换的过程，在转换过程中，需要注意 ID 的类型，以及一对多和多对多的情况。

1. #前面差异分析得到的res005文件
   
2. dta <- read.csv("res005.csv",header = T,row.names = 1,stringsAsFactors = F)
   
3. 
   
4. ## ------------------------------------------------------------------------
   
5. x <- rownames(dta)
   
6. length(x)
   
7. keytypes(org.Hs.eg.db)
   
8. eg <-  bitr(x, fromType="ENSEMBL", toType="ENTREZID", OrgDb="org.Hs.eg.db")
   
9. head(eg)
   
10. ids <- bitr(x, fromType="ENSEMBL", toType=c( "SYMBOL"), OrgDb="org.Hs.eg.db")
    
11. head(ids)
    
12. go <-  bitr(x,fromType = "ENSEMBL",toType = c("SYMBOL","GO","ONTOLOGY"),OrgDb = "org.Hs.eg.db")
    
13. head(go)

复制代码

三、GO 功能注释与富集
       Gene Ontology（简称 GO）是一个国际标准化的基因功能分类体系，提供了一套动态更新的标准词汇表（controlled vocabulary）来全面描述生物体中基因和基因产物的属性。GO 总共有三个 ontology（本体），分别描述基因的分子功能（molecular function）、细胞组分（cellular component）、参与的生物过程（biological process）。GO 的基本单位是 term（词条、节点），每个term 都对应一个属性。
基因 GO/KEGG 富集分析

1. GO 官网：

复制代码

1. gene <- eg$ENTREZID
   
2. gene.df <- bitr(gene, fromType = "ENTREZID",
   
3.                 toType = c("ENSEMBL", "SYMBOL"),
   
4.                 OrgDb = org.Hs.eg.db)
   
5. head(gene.df)
   
6. ggo <- groupGO(gene     = gene,
   
7.                OrgDb    = org.Hs.eg.db,
   
8.                ont      = "MF",
   
9.                level    = 3,
   
10.                readable = TRUE)
    
11. head(ggo)
    
12. View(as.data.frame(ggo))
    
13. 
    
14. ego <- enrichGO(gene          = gene,
    
15.                 OrgDb         = org.Hs.eg.db,
    
16.                 ont           = "CC",
    
17.                 pAdjustMethod = "BH",
    
18.                 readable      = TRUE)
    
19. head(ego)
    
20. View(as.data.frame(ego))

复制代码

四、GO功能注释可视化

1. barplot(ggo, drop=TRUE, showCategory=12)
   
2. 
   
3. ## ----fig.height=5, fig.width=8-------------------------------------------
   
4. barplot(ego, showCategory=15)
   
5. 
   
6. dotplot(ego)
   
7. 
   
8. goplot(ego)

复制代码

五、KEGG 功能注释与富集
       KEGG（Kyoto Encyclopedia of Genes and Genomes，京都基因与基因组百科全书）是一个关于基因组、酶促途径以及生物化学物质的在线数据库。虽然决定生物体基因分类的基因组测序工程有了飞速的发展，但对单个基因功能的研究仍然相差甚远。同时活细胞的生物学功能是许多分子相互作用的结果，不能仅仅归功于单个基因或单个分子。
KEGG 官网：http://www.genome.jp/kegg/

1. search_kegg_organism('ece', by='kegg_code')
   
2. search_kegg_organism('hsa', by='kegg_code')
   
3. ecoli <- search_kegg_organism('Escherichia coli', by='scientific_name')
   
4. dim(ecoli)
   
5. head(ecoli)
   
6. 
   
7. kk <- enrichKEGG(gene         = gene,
   
8.                  organism     = 'hsa',
   
9.                  pvalueCutoff = 0.05)
   
10. View(kk@result)
    
11. browseKEGG(kk, 'hsa04750')

复制代码

六、KEGG 功能注释可视化

1. barplot(kk)
   
2. dotplot(kk)
   
3. #自定义画图
   
4. y <- kk@result #保留kk结果
   
5. head(y)
   
6. par('mar')
   
7. par(mar=c(5.1,15.1,4.1,2.1))
   
8. barplot(height = y$Count[1:10],names.arg = y$Description[1:10],horiz = T,las=2)

复制代码

七、单独绘制 kegg 条目图

1. library(ggplot2)
   
2. pathway <-  read.csv("kegg.csv",header=T)
   
3. head(pathway)
   
4. colnames(pathway)
   
5. 
   
6. p <-  ggplot(data=pathway,mapping = aes(x=richFactor,y=Pathway))
   
7. p
   
8. p + geom_point()
   
9. p + geom_point(aes(size=AvsB))
   
10. p + geom_point(aes(size=AvsB,color=Qvalue))
    
11. p + geom_point(aes(size=AvsB,color=Qvalue)) + scale_colour_gradient(low="green",high="red")
    
12. p + geom_point(aes(size=AvsB,color=Qvalue)) + scale_colour_gradient(low="green",high="red")+labs(title="Top20 of pathway enrichment",x="Rich factor",y="Pathway Name",color="-log10(Qvalue)",size="Gene Numbers")
    
13. p + geom_point(aes(size=AvsB,color=Qvalue)) + scale_colour_gradient(low="green",high="red")+labs(title="Top20 of pathway enrichment",x="Rich factor",y="Pathway Name",color="-log10(Qvalue)",size="Gene Numbers")+theme_bw()

复制代码

      
       示例结果的ggplot2可视化