使用clusterProfiler進行GO富集分析
原創weixin_43569478 最后發布於2018-11-05 10:06:43 閱讀數 12741 收藏
展開
歡迎關注微信公眾號《生信修煉手冊》!
clusterProfiler是一個功能強大的R包，同時支持GO和KEGG的富集分析，而且可視化功能非常的優秀，本章主要介紹利用這個R包來進行Gene Ontology的富集分析。
進行GO分析時，需要考慮的一個基礎因素就是基因的GO注釋信息從何處獲取。Bioconductor上提供了以下19個物種的Org類型的包，包含了這些物種的GO注釋信息
packages organism
org.Ag.eg.db Anopheles
org.At.tair.db Arabidopsis
org.Bt.eg.db Bovine
org.Ce.eg.db Worm
org.Cf.eg.db Canine
org.Dm.eg.db Fly
org.Dr.eg.db Zebrafish
org.EcK12.eg.db E coli strain K12
org.EcSakai.eg.db E coli strain Sakai
org.Gg.eg.db Chicken
org.Hs.eg.db Human
org.Mm.eg.db Mouse
org.Mmu.eg.db Rhesus
org.Pf.plasmo.db Malaria
org.Pt.eg.db Chimp
org.Rn.eg.db Rat
org.Sc.sgd.db Yeast
org.Ss.eg.db Pig
org.Xl.eg.db Xenopus
clusterProfiler在做GO分析時，會讀取這些包中的GO信息，從這里也可以看出，Bioconductor的生態體系，不同R包的數據傳遞和依賴組成了其生態環境。
對於以上19個物種，只需要安裝對應的org包，clusterProfile就會自動從中獲取GO注釋信息，我們只需要差異基因的列表就可以了，使用起來非常方便。
在進行富集分析時，支持以下兩種算法
1. Over-Representation Analysis
過表達分析其實就是費舍爾精確檢驗，對於以上19個物種，分析的代碼如下
ego <- enrichGO(
gene = gene,
keyType = "ENTREZID",
OrgDb = org.Hs.eg.db,
ont = "CC",
pAdjustMethod = "BH",
pvalueCutoff = 0.01,
qvalueCutoff = 0.05,
readable = TRUE)
1
2
3
4
5
6
7
8
9
gene就是差異基因對應的向量，keyType指定基因ID的類型，默認為ENTREZID, 該參數的取值可以參考keytypes(org.Hs.eg.db)的結果， 建議采用ENTREZID, OrgDb指定該物種對應的org包的名字，ont代表GO的3大類別，BP, CC, MF; pAdjustMethod指定多重假設檢驗矯正的方法，cufoff指定對應的閾值，readable=TRUE代表將基因ID轉換為gene symbol。
2. Gene Set Enrichment Analysis
ego <- gseGO(
geneList = geneList,
OrgDb = org.Hs.eg.db,
ont = "CC",
nPerm = 1000,
minGSSize = 100,
maxGSSize = 500,
pvalueCutoff = 0.05,
verbose = FALSE)
1
2
3
4
5
6
7
8
9
GSEA分析通過置換檢驗來計算p值，nPerm指定置換次數。
對於這19個物種之外的其他物種，也支持讀取基因的GO注釋文件，然后進行分析，注釋文件的格式如下
GeneId GO Description
1 GO:0005819 spindle
2 GO:0072686 mitotic spindle
3 GO:0000776 kinetochore
只需要3列信息即可，第一列為geneID, 第二列為基因對應的GO編號，第三列為GO的描述信息。這3列的順序是無所謂的， 只要包含這3種信息就可以了。
讀取該文件，進行分析的代碼如下
data <- read.table(
"go_annotation.txt",
header = T,
sep = "\t")
go2gene <- data[, c(2, 1)]
go2name <- data[, c(2, 3)]
# 費舍爾精確檢驗
x <- enricher(
gene,
TERM2GENE = go2gene,
TERM2NAME = go2name)
# GSEA富集分析
x <- GSEA(
gene,
TERM2GENE = go2gene,
TERM2NAME = go2name)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
對於GO富集分析的結果，clusterProfiler提供了以下幾種可視化策略
1. barplot
用散點圖展示富集到的GO terms,用法如下
barplot(ego, showCategory = 10)
1
生成的圖片如下
橫軸為該GO term下的差異基因個數，縱軸為富集到的GO Terms的描述信息， showCategory指定展示的GO Terms的個數，默認展示顯著富集的top10個，即p.adjust最小的10個。注意的顏色對應p.adjust值，從小到大，對應藍色到紅色。
2. dotplot
用散點圖展示富集到的GO terms,用法如下
dotplot(ego, showCategory = 10)
1
生成的圖片如下
橫軸為GeneRatio, 代表該GO term下的差異基因個數占差異基因總數的比例，縱軸為富集到的GO Terms的描述信息， showCategory指定展示的GO Terms的個數，默認展示顯著富集的top10個，即p.adjust最小的10個。圖中點的顏色對應p.adjust的值，從小到大，對應藍色到紅色，大小對應該GO terms下的差異基因個數，個數越多，點越大。
3. GO有向無環圖
調用topGO來實現GO有向無環圖的繪制，代碼如下
plotGOgraph(ego)
1
生成的圖片如下
矩形代表富集到的top10個GO terms, 顏色從黃色過濾到紅色，對應p值從大到小。
4. goplot
和plotGOgraph展示的信息一樣，都是GO有向無環圖，只不過采用了igraph中的布局方式，用法如下
goplot(ego)
1
生成的圖片如下
5. emapplot
對於富集到的GO terms之間的基因重疊關系進行展示，如果兩個GO terms系的差異基因存在重疊，說明這兩個節點存在overlap關系，在圖中用線條連接起來，用法如下
emapplot(ego, showCategory = 30)
1
生成的圖片如下
每個節點是一個富集到的GO term, 默認畫top30個富集到的GO terms, 節點大小對應該GO terms下富集到的差異基因個數，節點的顏色對應p.adjust的值，從小到大，對應藍色到紅色。
6. cnetplot
對於基因和富集的GO terms之間的對應關系進行展示，如果一個基因位於一個GO Terms下，則將該基因與GO連線，用法如下
cnetplot(ego, showCategory = 5)
1
生成的圖片如下
圖中灰色的點代表基因，黃色的點代表富集到的GO terms, 默認畫top5富集到的GO terms, GO 節點的大小對應富集到的基因個數。更多用法和細節請參考官方文檔。
————————————————
版權聲明：本文為CSDN博主「weixin_43569478」的原創文章，遵循 CC 4.0 BY-SA 版權協議，轉載請附上原文出處鏈接及本聲明。
原文鏈接：https://blog.csdn.net/weixin_43569478/article/details/83744242