应用NMF替代层级聚类算法

日期：2021-02-04 浏览：

生信小白团留言板留言参加探讨，或是关心同名的手机微信微信公众号biotrainee

前边大家在实例教程：应用R包deconstructSigs依据己知的signature开展占比推论，圆满的把50八个患者，依据1一个signature开展了占比推论，获得的占比引流矩阵以一般的热点图，及其pheatmap包内置的层级聚类算法以下：

编码是：

rm(list=ls())
options(stringsAsFactors = F) 
load(file = 'mut.wt_from_denovo.Rdata')
a2=mut.wt
## 非负引流矩阵溶解鉴别驱动器性signature并制图
nmf.input - t(a2)
nmf.input - nmf.input[setdiff(rownames(nmf.input),"unknown"),] #除掉unknown
# 必须除去这些沒有被测算到的signatures
nmf.input - nmf.input[rowSums(nmf.input) 0,]
library(pheatmap)
pheatmap(nmf.input)

<_denovo_results.Rdata文档，来源于于实例教程：应用R包SomaticSignatures开展denovo的signature推论。

能看到，一部分患者的S1,S2占比较高，很显著，她们应当是一个组，随后剩下的患者能够依据S10开展高矮排序，其他signature仿佛是酱油。

那样的觉得，实际上便可令其用NMF优化算法来完成，特别是在是层级聚类算法其实不能非常好的把样版开展“渭泾分明”的排序。

第一步：分辨最好排序

便是数次运作 nmf ic系数：

library(NMF)
ranks - 2:10
estim - lapply(ranks, function(r){
 fit - nmf(nmf.input, r, nrun = 5, seed = 4, method = "lee") # nrun设定为5以防运作時间太长
 list(fit = fit, consensus = consensus(fit), .opt = "vp",coph = cophcor(fit))
names(estim) - paste('rank', ranks)
sapply(estim, '[[', 'coph')

png("ic coefficient for seleting optimal nmf rank.png")
par(cex.axis=1.5)
plot(ranks, sapply(estim, '[[', 'coph'), 
 xlab="", ylab="", type="b", 
 col="red", lwd=4,xaxt="n")
axis(side = 1, at=1:10)
title(xlab="number of clusters", ylab="ic coefficient", cex.lab=1.5)
dev.off()

以下所显示：

ic coefficient for seleting optimal nmf rank

ic评分挑选rank=4，由于它是转折点，可是零元，10钟头课堂教学视頻直播间《跟随百度搜索百度李彦宏学习培训恶性肿瘤遗传基因组测序数据信息剖析》这一参考文献里边设定为3。这儿，大家就重视参考文献結果吧。

第二步：挑选signature

前边大家的50八个患者，全是1一个signature，可是呢，大家的NMF优化算法运作之后，能看到有一些signature实际上对样版聚类算法排序并沒有实际意义，因此大家必须选择一下，编码以下：

rank - 3
seed - 2019620
rownames(nmf.input) - gsub("Signature","Sig",rownames(nmf.input)) # 行名简单化
mut.nmf - nmf(nmf.input, 
 rank = rank, 
 seed = seed, 
 method = "lee") 
index - extractFeatures(mut.nmf,"max") 
sig.order - unlist(index)

获得的結果以下；

 index
[[1]]
[1] 2 1
[[2]]
[1] 5 3 11 4
[[3]]
[1] 10 9
attr(,"method")
[1] "max"
 sig.order
[1] 2 1 5 3 11 4 10 9

能看到，nmf优化算法获得的3类，全是有各有侧重的signature，1一个signature也仅有八个必须事后开展剖析。以下图，能看到不一样nmf类有各有的侧重的signature。

basismap

第三步：应用选择出的signature再度NMF

 nmf.input2 - nmf.input[sig.order,]
library(pheatmap)
pheatmap(nmf.input2,cluster_rows = T,cluster_cols = F)
mut.nmf2 - nmf(nmf.input2, 
 rank = rank, 
 seed = seed, 
 method = "lee") 
group - predict(mut.nmf2) # 明确提出亚型
table(group)

这一情况下，患者的归类便是最终的nmf区别变成的3类。渭泾分明，以下：

consensusmap

番外：一些可视性化涵数

关键是再次参照每一个nmf类里边的不一样signature的占比，及其不一样nmf类的有关性热点图。

sample.order - names(group[order(group)])
#设定色调
jco - c("#2874C5","#EABF00","#C6524a","#868686")
png(file = "consensusmap.png")
consensusmap(mut.nmf2,
 labRow = NA,
 labCol = NA,
 annCol = data.frame("cluster"=group[colnames(nmf.input)]),
 annColors = list(cluster=c("1"=jco[1],"2"=jco[2],"3"=jco[3],"4"=jco[4])))
dev.off()
png(file = "basismap.png" ) 
# 从这幅图能够较为清楚地见到各亚型中的驱动器signature（深色），与下边的nmf heatmap相匹配
basismap(mut.nmf2,
 cexCol = 1,
 cexRow = 0.3,
 annColors=list(c("1"=jco[1],"2"=jco[2],"3"=jco[3],"4"=jco[4])))
dev.off()
aheatmap(as.matrix(nmf.input2[,sample.order]), 
 Rowv=NA, 
 Colv=NA, 
 annCol = data.frame("cluster"=group[sample.order]),
 annColors = list(cluster=c("1"=jco[1],"2"=jco[2],"3"=jco[3],"4"=jco[4])),
 color=c("#EAF0FA","#6081C3","#3454a7"), # 例文的深蓝色渐变色
 revC=TRUE, 
 cexCol = 0.3, 
 cexRow = 0.3,
 filename = "NMF_heatmap.pdf")
ac=data.frame(group=paste0('NMF_',group))
rownames(ac)=colnames(nmf.input2)
pheatmap(nmf.input2,show_colnames = F,annotation_col = ac)
pheatmap(nmf.input2[,sample.order],
 show_colnames = F,cluster_cols = F,annotation_col = ac)
pheatmap(nmf.input2,show_colnames = F,annotation_col = ac,
 filename = "NMF_heatmap1.pdf")
pheatmap(nmf.input2[,sample.order],
 show_colnames = F,cluster_cols = F,annotation_col = ac,
 filename = "NMF_heatmap2.pdf")
dev.off()

可视性化是一条就此止境的路。大伙儿给油哈！

网站知识

联系方式丨CONTACT

全国热线：18720358503
传真热线：18720358503
Q Q咨询：2639601583
企业邮箱：2639601583@qq.com

相关知识

应用NMF替代层级聚类算法

网站知识

联系方式丨CONTACT