今天繼續(xù)昨天推文的內(nèi)容，今天的內(nèi)容介紹如何在氣泡圖和左側(cè)和上方添加聚類樹(shù)圖，今天的內(nèi)容主要參考

https://mp.weixin.qq.com/s/XVl2MoOsT7pB1wNJmltoVw

這篇論文是在簡(jiǎn)書 土豆學(xué)生信 分享的內(nèi)容看到的刷后。簡(jiǎn)書的鏈接是 http://www.reibang.com/p/bbf9cb13b41a

論文是

image.png

論文對(duì)應(yīng)的代碼是公開(kāi)的 https://github.com/ajwilk/2020_Wilk_COVID

image.png

今天重復(fù)的內(nèi)容是論文中的figure2f

image.png

按照論文提供的代碼得到了畫圖用到的數(shù)據(jù)，部分?jǐn)?shù)據(jù)如下

image.png

但是用他提供的畫圖代碼沒(méi)有能夠畫出圖來(lái)贤姆。因?yàn)樗玫搅艘粋€(gè)dot_plot()函數(shù)煎源，沒(méi)有找到這個(gè)函數(shù)是怎么來(lái)的。既然已經(jīng)拿到了數(shù)據(jù)瘟忱，就用ggplot2自己來(lái)畫吧

讀入數(shù)據(jù)做氣泡圖奥额，

data.final<-read.csv("NM/figure2f.csv",header=T,check.names=F)
head(data.final)

library(ggplot2)
ggplot(data.final,aes(x=features.plot,y=id))+
  geom_point(aes(size=`Percent expressed`,
                 color=`Average expression`))+
  theme_bw()+
  theme(panel.grid = element_blank(),
        axis.text.x=element_text(angle=90,hjust = 1,vjust=0.5))+
  scale_color_gradient(low="lightgrey",high="blue")+
  labs(x=NULL,y=NULL)+
  guides(size=guide_legend(order=3))

image.png

以y軸為變量，做層次聚類访诱，并使用ggtree展示層次聚類結(jié)果

聚類用到的是平均表達(dá)量那一列

df<-data.final[,c(1,2,4)]

首先是長(zhǎng)格式數(shù)據(jù)轉(zhuǎn)換為寬格式

df1<-reshape2::dcast(df,id~features.plot,value.var = "Average expression")
rownames(df1)<-df1$id
df1.1<-df1[,2:22]

image.png

層次聚類垫挨，ggtree展示結(jié)果

df1.1.clust<-hclust(dist(df1.1))
df2.1.clust<-hclust(dist(df2.1))
library(ggtree)
p2<-ggtree(df1.1.clust)
p2+
  geom_tiplab()+
  xlim(NA,7)

image.png

使用aplot包拼圖

library(ggplot2)
p1<-ggplot(data.final,aes(x=features.plot,y=id))+
  geom_point(aes(size=`Percent expressed`,
                 color=`Average expression`))+
  theme_bw()+
  theme(panel.grid = element_blank(),
        axis.text.x=element_text(angle=90,hjust = 1,vjust=0.5))+
  scale_color_gradient(low="lightgrey",high="blue")+
  labs(x=NULL,y=NULL)+
  guides(size=guide_legend(order=3))

library(aplot)
p1%>%
  insert_left(p2,width = 0.2)

image.png

接下來(lái)就是在上方疊加聚類樹(shù)，一樣的操作

df2<-reshape2::dcast(df,features.plot~id,value.var = "Average expression")
rownames(df2)<-df2$features.plot
df2.1<-df2[,2:15]
df2.1.clust<-hclust(dist(df2.1))
p3<-ggtree(df2.1.clust)+
  #geom_tiplab(angle=90)+
  #theme_tree2()+
  layout_dendrogram()
p3
p1%>%
  insert_left(p2,width = 0.2)%>%
  insert_top(p3,height = 0.2)

這里多了一個(gè)知識(shí)點(diǎn)是ggtree作圖默認(rèn)開(kāi)口樹(shù)的方向是向右触菜，如果需要把開(kāi)口改成向下九榔，需要加上layout_dendrogram()函數(shù)

最終的結(jié)果如下

image.png

這里和論文中的圖有些不一致，可能是聚類算法的原因;ggtree有一個(gè)默認(rèn)的從上到下排序涡相，比如左側(cè)的樹(shù)現(xiàn)在第一個(gè)是H6,第二個(gè)是H5帚屉，如果想把H5放到第一個(gè)也是可以實(shí)現(xiàn)的，可以參考之前的推文 R語(yǔ)言ggtree按照指定的節(jié)點(diǎn)旋轉(zhuǎn)樹(shù)

歡迎大家關(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é)入門學(xué)習(xí)資料及自己的學(xué)習(xí)筆記！