多組差異基因富集分析結(jié)果聚類 | 獲得共表達(dá)term的富集圖

本教程問題:
1.mmGO.rdata文件的制作枣抱,是否不需要這一步,更改一下腳本即可實(shí)現(xiàn)呢蚊惯?
2.這樣的圖形美化愿卸,目前的圖還不是很滿意,需要進(jìn)行美化拣挪。
這些問題擦酌,需要我們大家一起解決,切分享哦菠劝!


一赊舶、前言

我們在做組學(xué)時,常常遇到多多個處理進(jìn)行差異分析赶诊,獲得差異基因或差異代謝物笼平。此后,是對這些差異結(jié)果進(jìn)行GO或KEGG富集分析舔痪。這是一個很普遍的結(jié)果寓调,但是這樣一弄后,我們對多個組合差異基因進(jìn)行富集分析后锄码,獲得多個GO或KEGG富集分析結(jié)果夺英,圖很多,但是重復(fù)的term也是很多滋捶,那如何闡述就是個問題痛悯,以及很多圖列在論文中也不是很美觀,只能間部分圖放在附件中重窟。



那么载萌,我們是不是可以間多個GO或KEGG的term進(jìn)行合并繪圖呢?這個想法在很久以前就像做,但是一直沒找到相關(guān)的教程扭仁。今天垮衷,突發(fā)奇想,無意間看到這樣的教程乖坠,那就重復(fù)一下吧搀突。

PS:在這個教程中,自己很多的參數(shù)還是沒有很好的了解瓤帚,如果你看過或是已經(jīng)可以深入了解描姚,希望你可以分享一下哦。


二戈次、圖形效果


這是最終圖形的效果,總體來說筒扒,還是基本滿足我們的需求怯邪,但是很多細(xì)節(jié)還是需要調(diào)整。

三花墩、開始繪圖

3.1 導(dǎo)入相關(guān)的包

## 導(dǎo)入R包
library(plyr)
library(stringr)
library(ape)
library(GOSemSim)
library(ggtree)   ## 進(jìn)化樹
library(scales)
library(cowplot)  ## 合圖
library(ggplot2)

設(shè)置路勁

setwd("D:\\小杜的生信筆記\\20220605GO富集聚類")

導(dǎo)入GO富集結(jié)果悬秉,這些結(jié)果在外面前期的教程中可獲得:clusterProfiler包 |GO、KEGG功能富集分析 | 值得收藏 冰蘑,GO和泌、KEGG功能富集分析 | 功能富集網(wǎng)絡(luò)圖、熱圖繪制(代碼重新)


多組富集分析結(jié)果的文件名統(tǒng)一命名祠肥,便于導(dǎo)入

fnames <- Sys.glob("enrichGO*.csv")
fnames

# GO term的合集
ego.ID <- unique(ego.all[,c(2:4)])
head(ego.ID)

> head(ego.ID)
           ID                                               Description   BgRatio  Bg
1  GO:0046777                               protein autophosphorylation 235/23239 235
2  GO:0034329                                    cell junction assembly 200/23239 200
3  GO:0034330                                cell junction organization 248/23239 248
5  GO:0098742 cell-cell adhesion via plasma-membrane adhesion molecules 199/23239 199
11 GO:0034332                            adherens junction organization 121/23239 121
12 GO:0043405                         regulation of MAP kinase activity 293/23239 293
```
## 3.3 GO term過濾
```
#刪掉超過300個基因的GO term
#如果不想刪除武氓,可以注釋掉下面這段
ego.ID$Bg <- as.numeric(str_split_fixed(ego.ID$BgRatio, "/",2)[,1]) #提取GO term包含的基因數(shù)量
ego.ID <- ego.ID[ego.ID$Bg < 300,] 
dim(ego.ID)

#刪掉少于100個基因的GO term
#具體怎樣篩選,要根據(jù)你自己的數(shù)據(jù)靈活設(shè)置
ego.ID <- ego.ID[ego.ID$Bg > 100,] 
dim(ego.ID)
head(ego.ID)
```
## 3.4 橫向合并多組富集結(jié)果仇箱,用于畫熱圖
```
## 橫向合并多組富集結(jié)果县恕,用于畫熱圖
MyMerge <- function(x, y){
  df <- merge(x, y, by= "ID", all.x= TRUE, all.y= TRUE)
  return(df)
}
ego.m <- Reduce(MyMerge, fdataset)
head(ego.m)
```
```
#只保留GO term ID和各組的p.adjust
ego.m <- ego.m[,c(1,4,7,10)] #此處有三組,如果有更多組剂桥,要在后面繼續(xù)加
#提取篩選過的GO term
ego.m <- merge(ego.ID[,1:2], ego.m, by= "ID", all.x= TRUE)
rownames(ego.m) <- ego.m$Description
ego.m$ID <- NULL
ego.m$Description <- NULL
```
```
#把列名改為組名
#colnames(ego.m)<- str_remove(fnames, ".csv")
colnames(ego.m) <- paste0("G", seq(1:length(fnames)))
head(ego.m)

> head(ego.m)
                                                       G1          G2          G3
tissue homeostasis                                     NA 0.008364684 0.004178739
regulation of cell-matrix adhesion                     NA 0.003506443 0.003381501
hematopoietic progenitor cell differentiation 0.007418708 0.002188681          NA
myeloid cell homeostasis                      0.012630541          NA          NA
myeloid leukocyte differentiation             0.003034455          NA          NA
regulation of leukocyte migration             0.005619600          NA 0.012970792

3.4 按GO term的相似性聚類

mmGO.rdata文件是這樣獲得的忠烛,但是我在考慮如果我們沒有自己.db是要如何獲得,是否可以根據(jù)我們前面的分析進(jìn)行獲得呢权逗?這個問題需要大家一起解決美尸,如果你知道可以分享一下,謝謝U遛薄师坎!*

library(org.Hs.eg.db)
hgGO <- godata('org.Hs.eg.db', ont="BP")
head(hgGO)
#save(hgGO, file="mmGO.rdata")

Q1,需要我們一起解決1伎选R倌汀!

#導(dǎo)入之前保存的文件
(load("mmGO.rdata")) 

#計算GO term之間的相似性
ego.sim <- mgoSim(ego.ID$ID, ego.ID$ID, semData=mmGO, measure="Wang", combine=NULL)
ego.sim[1:3, 1:3]

> ego.sim[1:3, 1:3]
                            protein autophosphorylation cell junction assembly
protein autophosphorylation                       1.000                  0.061
cell junction assembly                            0.061                  1.000
cell junction organization                        0.090                  0.683
                            cell junction organization
protein autophosphorylation                      0.090
cell junction assembly                           0.683
cell junction organization                       1.000

#用GO term作為行名、列名惶岭,便于查看和畫圖
rownames(ego.sim) <- ego.ID$Description
colnames(ego.sim) <- ego.ID$Description
ego.sim[1:3, 1:3]

繪制聚類圖

#給GO term分類
tree <- nj(as.dist(1-ego.sim))
ggtree(tree) + geom_tiplab() + #寫GO term
  geom_text2(aes(subset=!isTip, label=node), hjust=-.3) + #寫node編號
  coord_cartesian(xlim=c(-.1,1.3)) #左右兩側(cè)留出合適的空間

此圖節(jié)點(diǎn)就是我們后續(xù)分類的標(biāo)準(zhǔn)



3.5 繪圖

繪圖代碼在下面的教程中寿弱,如果獲得也可以在知乎、公眾號中進(jìn)行獲得代碼數(shù)據(jù)事例數(shù)據(jù)按灶。
多組差異基因富集分析結(jié)果聚類 | 獲得共表達(dá)term的富集圖


@小杜的生信筆記症革,主要發(fā)表或收錄生物信息學(xué)的教程,以及基于R的分析和可視化(包括數(shù)據(jù)分析鸯旁,圖形繪制等)噪矛;分享感興趣的文獻(xiàn)和學(xué)習(xí)資料!

?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
  • 序言:七十年代末,一起剝皮案震驚了整個濱河市铺罢,隨后出現(xiàn)的幾起案子艇挨,更是在濱河造成了極大的恐慌,老刑警劉巖韭赘,帶你破解...
    沈念sama閱讀 218,941評論 6 508
  • 序言:濱河連續(xù)發(fā)生了三起死亡事件缩滨,死亡現(xiàn)場離奇詭異,居然都是意外死亡泉瞻,警方通過查閱死者的電腦和手機(jī)脉漏,發(fā)現(xiàn)死者居然都...
    沈念sama閱讀 93,397評論 3 395
  • 文/潘曉璐 我一進(jìn)店門,熙熙樓的掌柜王于貴愁眉苦臉地迎上來袖牙,“玉大人侧巨,你說我怎么就攤上這事”薮铮” “怎么了司忱?”我有些...
    開封第一講書人閱讀 165,345評論 0 356
  • 文/不壞的土叔 我叫張陵,是天一觀的道長碉怔。 經(jīng)常有香客問我烘贴,道長,這世上最難降的妖魔是什么撮胧? 我笑而不...
    開封第一講書人閱讀 58,851評論 1 295
  • 正文 為了忘掉前任桨踪,我火速辦了婚禮,結(jié)果婚禮上芹啥,老公的妹妹穿的比我還像新娘锻离。我一直安慰自己,他們只是感情好墓怀,可當(dāng)我...
    茶點(diǎn)故事閱讀 67,868評論 6 392
  • 文/花漫 我一把揭開白布汽纠。 她就那樣靜靜地躺著,像睡著了一般傀履。 火紅的嫁衣襯著肌膚如雪虱朵。 梳的紋絲不亂的頭發(fā)上,一...
    開封第一講書人閱讀 51,688評論 1 305
  • 那天,我揣著相機(jī)與錄音碴犬,去河邊找鬼絮宁。 笑死,一個胖子當(dāng)著我的面吹牛服协,可吹牛的內(nèi)容都是我干的绍昂。 我是一名探鬼主播,決...
    沈念sama閱讀 40,414評論 3 418
  • 文/蒼蘭香墨 我猛地睜開眼偿荷,長吁一口氣:“原來是場噩夢啊……” “哼窘游!你這毒婦竟也來了?” 一聲冷哼從身側(cè)響起跳纳,我...
    開封第一講書人閱讀 39,319評論 0 276
  • 序言:老撾萬榮一對情侶失蹤忍饰,失蹤者是張志新(化名)和其女友劉穎,沒想到半個月后棒旗,有當(dāng)?shù)厝嗽跇淞掷锇l(fā)現(xiàn)了一具尸體喘批,經(jīng)...
    沈念sama閱讀 45,775評論 1 315
  • 正文 獨(dú)居荒郊野嶺守林人離奇死亡,尸身上長有42處帶血的膿包…… 初始之章·張勛 以下內(nèi)容為張勛視角 年9月15日...
    茶點(diǎn)故事閱讀 37,945評論 3 336
  • 正文 我和宋清朗相戀三年铣揉,在試婚紗的時候發(fā)現(xiàn)自己被綠了。 大學(xué)時的朋友給我發(fā)了我未婚夫和他白月光在一起吃飯的照片餐曹。...
    茶點(diǎn)故事閱讀 40,096評論 1 350
  • 序言:一個原本活蹦亂跳的男人離奇死亡逛拱,死狀恐怖,靈堂內(nèi)的尸體忽然破棺而出台猴,到底是詐尸還是另有隱情朽合,我是刑警寧澤,帶...
    沈念sama閱讀 35,789評論 5 346
  • 正文 年R本政府宣布饱狂,位于F島的核電站曹步,受9級特大地震影響,放射性物質(zhì)發(fā)生泄漏休讳。R本人自食惡果不足惜讲婚,卻給世界環(huán)境...
    茶點(diǎn)故事閱讀 41,437評論 3 331
  • 文/蒙蒙 一、第九天 我趴在偏房一處隱蔽的房頂上張望俊柔。 院中可真熱鬧筹麸,春花似錦、人聲如沸雏婶。這莊子的主人今日做“春日...
    開封第一講書人閱讀 31,993評論 0 22
  • 文/蒼蘭香墨 我抬頭看了看天上的太陽留晚。三九已至酵紫,卻和暖如春,著一層夾襖步出監(jiān)牢的瞬間,已是汗流浹背奖地。 一陣腳步聲響...
    開封第一講書人閱讀 33,107評論 1 271
  • 我被黑心中介騙來泰國打工橄唬, 沒想到剛下飛機(jī)就差點(diǎn)兒被人妖公主榨干…… 1. 我叫王不留,地道東北人鹉动。 一個月前我還...
    沈念sama閱讀 48,308評論 3 372
  • 正文 我出身青樓轧坎,卻偏偏與公主長得像,于是被迫代替她去往敵國和親泽示。 傳聞我的和親對象是個殘疾皇子缸血,可洞房花燭夜當(dāng)晚...
    茶點(diǎn)故事閱讀 45,037評論 2 355

推薦閱讀更多精彩內(nèi)容