跟着小果Get无监督聚类!-ConsensusClusterPlus

小伙伴在做生信分析的时候,会经常见到聚类这个名词,却不知道是什么意义。其实平时我们做基因热图时候就用到了聚类的思想。有的小伙伴在做癌症分析时候,也见过给疾病做分子亚型,这个时候无监督聚类就派上用场了。

所以今天小果带大家了解无监督聚类,并学习如何使用他。

首先我们先大致了解一下什么是无监督聚类。就好比一致性聚类方法包括从一组项目中进行次抽样,例如微阵列,并确定特定簇计数(k)的簇。然后,对共识值,两个项目占在同一子样本中发生的次数中有相同的聚类,计算并存储在每个k的对称一致矩阵中。

具体的含义我们不在此一一赘述。

小果在这里用到的是ConsensusClusterPlus这个包对基因进行聚类。

所需要的是下面这种格式的基因表达谱:

#install.packages(“BiocManager”)

#BiocManager::install(‘ConsensusClusterPlus’)#小伙伴没有的自己要下载一下

library(tidyverse)

library(ConsensusClusterPlus)

#导入数据

exp <- read.table(“text1.txt”,sep = “\t”,row.names = 1,check.names = F,stringsAsFactors = F,header = T)

d=as.matrix(exp)

gene <- c(“MMP9″,”CCND1″,”STAT1″,”AR”,”AKR1C3″,”GSTP1″,”PGR”)

d <- d[gene,]

mads=apply(d,1,mad)

d=d[rev(order(mads))[1:7],] #7是基因个数

d = sweep(d,1, apply(d,1,median,na.rm=T))

title=(“JULEI”) ##文件夹输出图片的位置

set.seed(1) #我发现设不设置种子都一样

#这里小果选取自己的示例数据,小伙伴们整理自己想去聚类的基因。

然后我们打开JULEI这个文件夹会发现得到了聚类的结果:
一致性矩阵

一致性累积分布函数图

碎石图

在这里小果就不全部展示了,小伙伴们自己做的话可以慢慢看。

从上述碎石图看出呢,在k=7,k=8左右时候呢,比较平稳,前面K=3,和k=4这两个分组数间相对差值变化越大,说明该分组越有价值。小伙伴需要对症下药,查看自己聚类的结果进行解释和后续的分析。

图中K值呢就代表簇,

接下来我们绘制另一组图:

聚类一致性(cluster-consensus)以及样品一致性(item-consensus)

title=title,clusterAlg=”hc”,distance=”pearson”,seed=1,plot=”pdf”)

results[[2]][[“consensusMatrix”]][1:5,1:5]

results[[2]][[“consensusTree”]]

results[[2]][[“consensusClass”]][1:5]

icl = calcICL(results,title=title,plot=”pdf”) ##画另一组图片

结果如下:

聚类一致性(cluster-consensus)

样品一致性(item-consensus)

小伙伴们把需要那些图片进行展示就可以了。

#上述就是绘制另一组图,然而我们还需要做的是把结果保存下来。

group<-results[[2]][[“consensusClass”]]

group<-as.data.frame(group)

group$group <- factor(group$group,levels=c(1,2))

save(group,file = “group_AY.Rda”)#这里把分组结果保存

然而,又会有小伙伴会问如何展示聚类后的效果呢。

小果在这里教大家一种方式用热图的方式进展展示,

绘制ConsensusClusterPlus后的热图:

这个时候就用到上述咱们保存的分组信息了:

跟着小果的代码往下跑就可以啦!

library(pheatmap)#调用一下热图的包

group <- group %>%

rownames_to_column(“sample”)

annotation <- group %>% arrange(group) %>% column_to_rownames(“sample”)

a <- group %>% arrange(group) %>% mutate(sample=substring(.$sample,1,12))

b <- t(exp_gene) %>%

as.data.frame() %>%

rownames_to_column(“sample”) %>%

mutate(sample=substring(.$sample,1,12))

c <- inner_join(a,b,”sample”) %>% .[,-2] %>% column_to_rownames(“sample”) %>% t(.)

pheatmap(c,annotation = annotation,

cluster_cols = F,fontsize=5,fontsize_row=5,

scale=”row”,show_colnames=F,

fontsize_col=3)

bk <- c(seq(-5,-0.1,by=0.01),seq(0,5,by=0.01))

pheatmap(c,annotation = annotation,

annotation_colors = list(group = c(“1″ =”#8B1C62”,”2″= “#EE9572″)),

cluster_cols = F,fontsize=10,fontsize_row=10,

scale=”row”,show_colnames=F,cluster_row = F,

color = c(colorRampPalette(colors = c(“#551A8B”,”white”))(length(bk)/2),colorRampPalette(colors = c(“white”,”#CD3700″))(length(bk)/2)),

#legend_breaks=seq(-3,3,1),

fontsize_col=3)

dev.off()

这里热图的参数比如:颜色的绘制都可以根据小伙伴自己的需求来修改,

效果就出来了,可以看出来在MMP9和CCND1等这几个基因聚类效果比较明显,小果在这里展示的是K=2情况下聚类的热图。小伙伴们自己选择K的值进行绘制。

下面就可以根据聚类的结果去后续的分析了。比如;分组的ssGSEA等等

小伙伴们需要多多学习代码的含义,以及对图片的理解才能,更好的去运用。

上述就是对无监督聚类法的学习。小伙伴们是否学会了呢