GOplot那個(gè)R語(yǔ)言包畫(huà)弦圖展示GO富集分析的結(jié)果吨灭,圖還挺漂亮的,但是如果是自己的數(shù)據(jù)到最終他的畫(huà)圖需要的輸入數(shù)據(jù)還挺麻煩的刑巧,我這邊寫(xiě)了一個(gè)python腳本喧兄,希望可以簡(jiǎn)化畫(huà)圖數(shù)據(jù)的準(zhǔn)備過(guò)程,腳本寫(xiě)的比較丑啊楚,但是湊合著能用
你自己需要準(zhǔn)備的數(shù)據(jù)是
- GO富集分析的結(jié)果
- 感興趣的基因列表(每行一個(gè)基因名)
- 感興趣的基因?qū)?yīng)的數(shù)據(jù)(比如表達(dá)量或者其他)
一個(gè)基因名對(duì)應(yīng)著一個(gè)數(shù)據(jù)吠冤,每行一個(gè) - 感興趣的Term的名字(每行一個(gè))
運(yùn)行腳本
python .\prepare_input_df_for_R_GOplot_GOchord.py .\genes.txt .\process.txt .\GO_enrich_result.tsv a b 3 4 output-1.txt
- 第一個(gè)位置是感興趣的基因名文件
- 第二個(gè)位置是感興趣的Term的文件
- 第三個(gè)位置是GO富集分析的結(jié)果文件
- 第四個(gè)位置是GO富集分析結(jié)果文件的分隔符
a代表制表符
b代表逗號(hào)
d代表斜線/ - 第五個(gè)位置是GO富集分析結(jié)果文件中基因名那一列的分隔符 abd同上
- 第六個(gè)是GO富集分析結(jié)果文件中Term對(duì)應(yīng)的是哪一列
- 第七個(gè)是GO富集分析結(jié)果文件中基因名對(duì)應(yīng)的是哪一列
- 第八個(gè)參數(shù)是輸出文件的名字
這樣就獲得了output-1.txt
這個(gè)文件,接下來(lái)是R語(yǔ)言里的代碼
df<-read.csv("D:/Jupyter/GOplot/example/output-1.txt",
header = T,
sep="\t",
row.names = 1,
check.names = F)
head(df)
df<-df[,1:dim(df)[2]-1]
首先是讀入數(shù)據(jù)恭理,python腳本輸出的文件是制表符分隔的
結(jié)尾處有一列空值拯辙,python里不知道如何刪掉,這里讀入數(shù)據(jù)以后再刪掉吧
which(rowSums(df) == 0)
數(shù)據(jù)每行都是0颜价,和1涯保,還得檢查一下有沒(méi)有哪一行全是0,如果全是0這一行得刪掉
df1<-read.csv("D:/Jupyter/GOplot/example/geneslogfc.csv")
df$logFC<-df1$logFC
讀取帶數(shù)值的基因文件周伦,這個(gè)基因文件的順序需要和感興趣的基因那個(gè)文件完全保持一致
最后就是畫(huà)圖了
library(GOplot)
library(ggplot2)
p1<-GOChord(as.matrix(df),
space = 0.02,
gene.order = 'logFC',
gene.space = 0.25, gene.size = 5)
ggsave(filename = "GOplot_example/p3.pdf",
p1,
width = 15,height = 15)
后記
python里是可以調(diào)用R語(yǔ)言的夕春,好好研究一些,看能不能直接一個(gè)python腳本出圖专挪,那樣就非常方便了
歡迎大家關(guān)注我的公眾號(hào)
小明的數(shù)據(jù)分析筆記本
小明的數(shù)據(jù)分析筆記本 公眾號(hào) 主要分享:1及志、R語(yǔ)言和python做數(shù)據(jù)分析和數(shù)據(jù)可視化的簡(jiǎn)單小例子片排;2、園藝植物相關(guān)轉(zhuǎn)錄組學(xué)速侈、基因組學(xué)率寡、群體遺傳學(xué)文獻(xiàn)閱讀筆記;3倚搬、生物信息學(xué)入門(mén)學(xué)習(xí)資料及自己的學(xué)習(xí)筆記冶共!