用SOM畫一個我認為還行的圖

用SOM package畫一個還行(我認為)的圖

SOM package 和kohonen package都可以做SOM cluster分析。kohonen的圖形展示要漂亮很多,還有ggsom包進一步優(yōu)化。我沒找到SOM包的圖形優(yōu)化包,決定自己畫一個。

SOM package

require(som)
SOMdata=data.frame(read.table("xxx.xls", header = TRUE))
data_som=SOMdata[,2:7] #這一步為了去掉rownames
#data_som.f <- filtering(as.matrix(data_som)) 這步
scaledata<-scale(data_som)
normaldata<-normalize(as.matrix(data_som),byrow = T)
data_som=som(normaldata,xdim=5, ydim=6, init="linear", alpha=NULL, 
                      alphaType="inverse", neigh="gaussian", topol="rect")
plot(data_som)

plot(data_som)這一步會得到一張中心數(shù)據(jù)點趨勢線圖,但是我不想跟大家都一樣,所以決定用ggplot2把數(shù)據(jù)重新畫一遍。中心數(shù)據(jù)在xx$code中

Make datasheet for ggplot2

得到的數(shù)據(jù)list(data_som)不能直接用ggplot2畫圖,必須轉(zhuǎn)化為dataframe

#make a dataframe for ggplot2
data_som_frame<- as.data.frame(data_som$code)
colnames(data_som_frame)<-c("A", "B", "C", "D", "E", "F")
data_som_frame$sub<-c(1:30)  #30是是SOM subclusters的數(shù)目= xdim × ydim
data_som_frame$nobs<-data_som$code.sum$nobs  #nobs是每個subcluster中gene的個數(shù)
#reshape datasheet into long-sheet that ggplot2 can handle
require(reshape2)
data_long<- melt(data_som_frame, id = c("sub","nobs"), variable.name = "condition", value.name = "code")

然后就可以用ggplot2想畫什么樣畫什么樣啦

p<-ggplot(data_long, aes(x=condition, y=code, group= sub))+geom_point(pch = 1,
         size = 2)+ scale_shape(solid=FALSE)+geom_line() 
p<-p+facet_wrap(~sub,labeller = "label_both")

labeller

其實到這一步就可以了,會優(yōu)化的繼續(xù)優(yōu)化。我不太會,我想把subcluster(1~30)標在圖上,同時顯示每個分組中gene的數(shù)目,意味著我需要在每一個facet的label上引用兩個單元格的內(nèi)容,另外也不想要原始的灰框,我采用了下面的方法:

p<-p+facet_wrap(~sub,labeller = "label_both")+theme(strip.background = element_blank())
p<-p+geom_text(aes(label = nobs), x = Inf, y = Inf, 
               hjust = 1.5, vjust = 3)
p

另外有一篇帖子 Put multi-variable facet_wrap labels on one line 這篇看起來非常簡潔,也可以。

p2<-ggplot(data_long, aes(x=condition, y=code, group=sub)) +
  geom_line()+ geom_point(size = 2) + facet_wrap(~paste(sub,nobs, sep="-"))
p2

最后用ggthemer 美化一下

require(ggthemes)
ggthemr('pale', layout = 'clean', spacing = 0.8, type = 'inner')
p
p2

最后是關于ggplot2labeller()的使用,有機會應該仔細研讀啊……畫個圖愁死了


更新:

用log2FoldChange后的數(shù)據(jù)畫圖

除了采用中心數(shù)據(jù)(xx$code),文獻中很多圖使用log2FoldChange來表示gene的變化,我也想用這種方法(而非code)來表示每個基因在每個條件下的變化,用一條平均數(shù)線來表示subcluster整體的變化,然后用hline()畫出零點水平直線。

#SOM-graph using original data (log2foldchange)
#將分組以indices的方式賦回原數(shù)據(jù)組,data222.log是原值經(jīng)過log2變換后的矩陣表
#data_som$visual$x[i]是SOM算法聚類后相應基因歸屬的組(x,y)
data222.log$x = c()
for(i in 1:2336) {data222.log$x[i] <-c(data_som$visual$x[i])}
data222.log$y = c()
for(i in 1:2336) {data222.log$y[i] <-c(data_som$visual$y[i])}

#reshape datasheet into long-sheet that ggplot2 can handle
data222.log.long<- melt(data222.log, id = c("x","y","gene"), 
                  variable.name = "condition", value.name = "code")

#ggplot2+hline()
p<-ggplot(data222.log.long, aes(x=condition, y=code, group=gene)) +
  geom_line(alpha = 3/5, colour = "#f3c57b")+geom_hline(aes(yintercept=0), colour="black") 

#加入平均值,加粗
p+ stat_summary(aes(group=x+y), fun.y=mean, geom="line", colour="#db735c",size=1)
    + facet_wrap(~paste(x,y, sep="-"),scales = "free")
p
最后編輯于
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成,瀏覽時請結(jié)合常識與多方信息審慎甄別。
平臺聲明:文章內(nèi)容(如有圖片或視頻亦包括在內(nèi))由作者上傳并發(fā)布,文章內(nèi)容僅代表作者本人觀點,簡書系信息發(fā)布平臺,僅提供信息存儲服務。

相關閱讀更多精彩內(nèi)容

  • 上周的暑期生信黑馬培訓有老師提出要做SOM分析,最后卡在code plot只能出segment plot卻出不來l...
    生信寶典閱讀 7,623評論 7 9
  • 工欲善其事,必先利其器??偨Y(jié)一下,方便多了。R語言還是很牛逼的,可以干很多事情。有一把順手的刀還是很重要的。 0....
    Liam_ml閱讀 4,861評論 1 60
  • Swift1> Swift和OC的區(qū)別1.1> Swift沒有地址/指針的概念1.2> 泛型1.3> 類型嚴謹 對...
    cosWriter閱讀 11,621評論 1 32
  • 今年剛剛畢業(yè),自己的專業(yè)是一名白衣天使,一直沒有出遠門的自己,想去外面的世界看看,所以選擇了遙遠的新疆。 ...
    夢芝_1e5d閱讀 373評論 0 1
  • 文/瑾初時 幾天前,一枚小仙女在周杰倫的演唱會上火了。 (演唱會現(xiàn)場的點歌環(huán)節(jié)): 杰倫:你要點什么歌曲? 小仙女...
    瑾初時_nunu閱讀 496評論 2 2

友情鏈接更多精彩內(nèi)容