clusterProfiler是一個功能強大的R包，同時支持GO和KEGG的富集分析，而且可視化功能非常的優(yōu)秀，本章主要介紹利用這個R包來進行Gene Ontology的富集分析。
進行GO分析時，需要考慮的一個基礎(chǔ)因素就是基因的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的生態(tài)體系，不同R包的數(shù)據(jù)傳遞和依賴組成了其生態(tài)環(huán)境。
對于以上19個物種，只需要安裝對應(yīng)的org包，clusterProfile就會自動從中獲取GO注釋信息，我們只需要差異基因的列表就可以了，使用起來非常方便。
在進行富集分析時，支持以下兩種算法
過表達分析其實就是費舍爾精確檢驗，對于以上19個物種，分析的代碼如下
ego <- enrichGO( gene = gene, keyType = "ENTREZID", OrgDb = org.Hs.eg.db, ont = "CC", pAdjustMethod = "BH", pvalueCutoff = 0.01, qvalueCutoff = 0.05, readable = TRUE)gene就是差異基因?qū)?yīng)的向量，keyType指定基因ID的類型，默認為ENTREZID, 該參數(shù)的取值可以參考keytypes(org.Hs.eg.db)的結(jié)果， 建議采用ENTREZID, OrgDb指定該物種對應(yīng)的org包的名字，ont代表GO的3大類別，BP, CC, MF; pAdjustMethod指定多重假設(shè)檢驗矯正的方法，cufoff指定對應(yīng)的閾值，readable=TRUE代表將基因ID轉(zhuǎn)換為gene symbol。
ego <- gseGO( geneList = geneList, OrgDb = org.Hs.eg.db, ont = "CC", nPerm = 1000, minGSSize = 100, maxGSSize = 500, pvalueCutoff = 0.05, verbose = FALSE)GSEA分析通過置換檢驗來計算p值，nPerm指定置換次數(shù)。
對于這19個物種之外的其他物種，也支持讀取基因的GO注釋文件，然后進行分析，注釋文件的格式如下
| GeneId | GO | Description |
|---|---|---|
| 1 | GO:0005819 | spindle |
| 2 | GO:0072686 | mitotic spindle |
| 3 | GO:0000776 | kinetochore |
只需要3列信息即可，第一列為geneID, 第二列為基因?qū)?yīng)的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)對于GO富集分析的結(jié)果，clusterProfiler提供了以下幾種可視化策略
用散點圖展示富集到的GO terms,用法如下
barplot(ego, showCategory = 10)生成的圖片如下
橫軸為該GO term下的差異基因個數(shù)，縱軸為富集到的GO Terms的描述信息， showCategory指定展示的GO Terms的個數(shù)，默認展示顯著富集的top10個，即p.adjust最小的10個。注意的顏色對應(yīng)p.adjust值，從小到大，對應(yīng)藍色到紅色。
用散點圖展示富集到的GO terms,用法如下
dotplot(ego, showCategory = 10)生成的圖片如下
橫軸為GeneRatio, 代表該GO term下的差異基因個數(shù)占差異基因總數(shù)的比例，縱軸為富集到的GO Terms的描述信息， showCategory指定展示的GO Terms的個數(shù)，默認展示顯著富集的top10個，即p.adjust最小的10個。
圖中點的顏色對應(yīng)p.adjust的值，從小到大，對應(yīng)藍色到紅色，大小對應(yīng)該GO terms下的差異基因個數(shù)，個數(shù)越多，點越大。
調(diào)用topGO來實現(xiàn)GO有向無環(huán)圖的繪制，代碼如下
plotGOgraph(ego)生成的圖片如下
和plotGOgraph展示的信息一樣，都是GO有向無環(huán)圖，只不過采用了igraph中的布局方式，用法如下
goplot(ego)生成的圖片如下
對于富集到的GO terms之間的基因重疊關(guān)系進行展示，如果兩個GO terms系的差異基因存在重疊，說明這兩個節(jié)點存在overlap關(guān)系，在圖中用線條連接起來，用法如下
emapplot(ego, showCategory = 30)生成的圖片如下
每個節(jié)點是一個富集到的GO term, 默認畫top30個富集到的GO terms, 節(jié)點大小對應(yīng)該GO terms下富集到的差異基因個數(shù)，節(jié)點的顏色對應(yīng)p.adjust的值，從小到大，對應(yīng)藍色到紅色。
對于基因和富集的GO terms之間的對應(yīng)關(guān)系進行展示，如果一個基因位于一個GO Terms下，則將該基因與GO連線，用法如下
cnetplot(ego, showCategory = 5)生成的圖片如下
圖中灰色的點代表基因，黃色的點代表富集到的GO terms, 默認畫top5富集到的GO terms, GO 節(jié)點的大小對應(yīng)富集到的基因個數(shù)。
更多用法和細節(jié)請參考官方文檔。
·end·
