R語言學(xué)習(xí) - 富集分析泡泡圖 剛一出品，Y叔就說有硬傷。Y叔是著名富集分析軟件clusterprofiler的原創(chuàng)，而且軟件內(nèi)集成dotplot, enrichmap,cnetmap (后續(xù)也實(shí)現(xiàn)這兩個(gè)的一步出圖)等畫圖方法，具體看這個(gè)教程 http://guangchuangyu.github.io/2016/01/go-analysis-using-clusterprofiler/ 或Biobabble公眾號(hào)。

這個(gè)意見得重視，不過大夏天的，鍋還是不能背著，回應(yīng)下Y叔的回應(yīng)。

一個(gè)command出圖，小白上點(diǎn)心可以的

圖是ggplot2畫的 (散點(diǎn)圖)，在腳本功能描述里有寫。只需要在終端輸入腳本的名字sp_enrichmentPlot.sh (腳本需在環(huán)境變量中) 就可以看到腳本功能描述、輸入文件樣例、參數(shù)列表和使用例子。

目前腳本的輸入是不支持Windows的Excel格式的，只支持TAB分割的文本文件。之前一直在Linux下工作，文件名都不寫后綴的。后來跟Windows用戶打交道多，為了方便他們打開，強(qiáng)行加了xls后綴，這一點(diǎn)有誤導(dǎo)。

腳本依賴R，這是應(yīng)該做一個(gè)判斷和提示的。但有R卻提示ggplot2包或其它包不存在，是可以用-i TRUE來自動(dòng)安裝包的。(最近在跟德春合作，完善包的自動(dòng)檢測(cè)和安裝，最后一起整合在分發(fā)包。)

Y叔的這兩點(diǎn)確實(shí)戳中了當(dāng)前文檔信息不完善的弱點(diǎn)，也希望大家多反映使用過程中的問題，幫助我們改進(jìn)文檔。

不過話說回來，如果仔細(xì)看了當(dāng)前腳本的功能描述和參數(shù)提示，上點(diǎn)心得小白不只可以一步畫圖，還可以隨意調(diào)整樣式。

不只clusterprofiler的用戶需要

如果用clusterprofiler做富集分析，write.table輸出結(jié)果，那么輸入文件、R、ggplot2都有了，一步出圖沒有問題。(Windows下沒bash？https://mp.weixin.qq.com/s/g1twNEsPWZb_tZFDyoTqVQ，看這里)

Y叔說既然有了clusterprofiler，一步步在R里面做不是很好，為啥要把數(shù)據(jù)導(dǎo)出來？

首先，轉(zhuǎn)換數(shù)據(jù)、存儲(chǔ)數(shù)據(jù)這一步是獨(dú)立于作圖的。圖只能顯示部分?jǐn)?shù)據(jù) (這點(diǎn)Y叔也有提到，用了simplify會(huì)好些，但也還會(huì)有不少通路)，所有富集條目導(dǎo)出作為文章附表，以顯示信息的全面和真實(shí)。

其次，即便可以畫出所有富集數(shù)據(jù) (用一副大圖)，也會(huì)先對(duì)結(jié)果做下篩選，一些特別基礎(chǔ)的、極父層的生物富集通路也會(huì)選擇不展示，優(yōu)先展示樣品屬性更相關(guān)的。畢竟是看在哪些通路里面富集，不是看在哪些通路里面最富集。所以需要導(dǎo)出數(shù)據(jù)，做下篩選，然后一步繪圖。

再次，不少人拿到的是已經(jīng)做好的富集分析結(jié)果 (可能是clusterprofiler的，也可能是其它軟件的)，總不能再導(dǎo)入clusterprofiler繪圖吧，可以有更簡(jiǎn)單的方式的。

最后，一個(gè)圖做出來，需要反復(fù)修改。今天用clusterprofiler做了富集分析，運(yùn)用dotplot配合不同參數(shù)出了圖；過幾天心情一變，想換個(gè)風(fēng)格了，怎么辦？再運(yùn)行一次clusterprofiler還是加載之前存儲(chǔ)的.Rdata。好像都不太方便，還是用導(dǎo)出的文本一步出圖吧。

關(guān)于硬傷

示例圖中沒有overlap不是繪圖腳本的問題，是數(shù)據(jù)篩選的問題。腳本如實(shí)反映篩選出的數(shù)據(jù)。如果選擇的通路在一個(gè)基因集合富集，另外一個(gè)不富集，就不顯示；如果選擇的通路是在多個(gè)基因集中都富集的，但富集程度不同，是可以顯示的; 而且這些通路會(huì)優(yōu)先顯示在圖的上部，下面我會(huì)給出例子。

一步出圖也可以定制

說到互動(dòng)，一步出圖不只可以，而且還可以記錄互動(dòng) (給每個(gè)輸出文件加個(gè)唯一的字符串做為文件名的一部分就可以了，不過這個(gè)沒用到，之前就沒寫這個(gè)參數(shù))。

拿Y叔的數(shù)據(jù)做個(gè)例子

library(clusterProfiler) library(dplyr) data(gcSample) x <- compareCluster(gcSample, 'enrichDO') # 獲取數(shù)據(jù)，取出每個(gè)組最富集的10個(gè)條目，存儲(chǔ)起來 # 不麻煩的吧 x <- x@compareClusterResult y = x %>% group_by(Cluster) %>% top_n(-10, pvalue) y = x[x$Description %in% y$Description,]

在樣品名中包含差異基因的數(shù)目

# 這個(gè)不解釋了，前面繪圖教程中提到過 y$Sample <- paste0(y$Cluster, "(",unlist(lapply(strsplit(y$GeneRatio, '/'), function(x) x[2])),")") # 這次存成tsv，省得誤會(huì)了 write.table(y, file="enrichment.tsv", sep="\t", col.names=T, row.names=F, quote=F)

首先出第一個(gè)圖，很簡(jiǎn)單，一個(gè)命令，設(shè)置下顏色。

+ 更正下，我這里作弊了一下，之前點(diǎn)的大小都是Gene Count這樣的純數(shù)字列來顯示的，不支持GeneRatio這樣的分?jǐn)?shù)形式的列 (上一版本，如果GeneRatio出現(xiàn)在橫軸是支持的)，做了下修改，也支持了。如果想用之前的代碼重復(fù)，把-s GeneRatio改為-s Count。 # 具體參數(shù)解釋看上一篇公眾號(hào)文章 # -C: 指定顏色，默認(rèn)是紅綠，可以是任何兩個(gè)R支持的顏色名字 sp_enrichmentPlot.sh -f enrichment.tsv -o Sample -T string -v Description -c p.adjust -s GeneRatio -C "'red','blue'" -w 16

圖略有不同，是因?yàn)榕判虻姆绞讲煌?。sp_enrichmentPlot.sh設(shè)置的是先按照出現(xiàn)在不同樣品最多的條目?jī)?yōu)先的策略，可以清晰的看到哪些是不同基因集共同富集的，哪些是不同基因集特異富集的。其它的都一樣。當(dāng)然沒有設(shè)置字體參數(shù)，字體略小了些。

說還想要個(gè)紫綠的，也可以，只需改兩個(gè)單詞，而不是加一段代碼scale_color_continuous(low='purple', high='green')。

sp_enrichmentPlot.sh -f enrichment.tsv -o Sample -T string -v Description -c p.adjust -s GeneRatio -C "'purple','green'" -w 16

改變下形狀，也可以啊。程序中都留有一個(gè)參數(shù)可以寫入ggplot2語句，就是下面的-z。那么，ggplot2能修改的樣式也都可以。常用的修改會(huì)做成一個(gè)參數(shù)，不常用的就只能直接寫命令了。能在Rstudio里面敲，那么也能在命令行敲。

#-z: 后跟任何合法的ggplot2語句 sp_enrichmentPlot.sh -f enrichment.tsv -o Sample -T string -v Description -c p.adjust -s GeneRatio -C "'purple','green'" -w 16 -z "+aes(shape=GeneRatio>0.1)"

-z “+ any legal ggplot2 command”

當(dāng)然這一步也可以不通過增加ggplot2語句來實(shí)現(xiàn)，在數(shù)據(jù)中加一列就好了，通過-S指定參數(shù)。

不加ggplot2語句, 只修改命令行參數(shù)，sp_enrichmentPlot.sh也可以調(diào)整圖例的位置、輸出文件的格式、是否對(duì)p.adj取log、是否分面等。

一步作圖的優(yōu)勢(shì)

一步作圖相比于直接寫R代碼有什么好處呢？

模塊化好。也就是Y叔提到的數(shù)據(jù)處理和可視化分開；一步作圖，只是作圖，不做無關(guān)的處理，更隨意。

易用性強(qiáng)。敲代碼，總不如給改數(shù)來得快。

重用性強(qiáng)。假如我要做十幾個(gè)分開的基因集合呢？一段段復(fù)制代碼？改錯(cuò)了或忘記改某個(gè)地方了怎么辦？

快速出圖。先快速出個(gè)原型，再接著調(diào)整。

生信寶典，一起換個(gè)角度學(xué)生信

http://mp.weixin.qq.com/s?__biz=MzI5MTcwNjA4NQ==&mid=2247484063&idx=1&sn=f4e93d428e4910b4abbee9c0430cd170&chksm=ec0dc715db7a4e0318b388ba2ab3d51677741421c42ada474a0ac6046a0699283014eae84b6f#rd

總結(jié)

以上是生活随笔為你收集整理的R 学习 - 富集分析泡泡图回应的全部?jī)?nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。