基于多组学数据的肿瘤亚型分组一直是研究者关注的一个方向,今天给大家介绍一个基于多组学进行聚类的R包SNFtool。此工具主要是利用相似网络融合将一个网络的多个视图融合在一起,构建一个整体的状态矩阵。算法的输入可以是特征向量、成对距离或成对相似度。学习到的状态矩阵可以用于检索、聚类和分类。其中核心的算法是谱聚类算法,其主要是完成对图的分割,找到最好的分割方式,来将图分割开来。这些图在这里就是我们构建的距离矩阵。简而言之,谱聚类先降维(特征分解),然后在低维空间用其它聚类算法(如KMeans、模糊聚类)进行聚类。
首先我们看下包的按装:
代码语言:javascript复制install.packages(“SNFtool”)
接下来我们利用包自带的数据进行整个流程的再现:
代码语言:javascript复制##数据的载入
data("Data1")
data("Data2")
data("dataL")
代码语言:javascript复制##数据的标准化,归一化为均值为0,标准差为1的矩阵
Data1 = standardNormalization(Data1);
Data2 = standardNormalization(Data2);
代码语言:javascript复制##欧氏距离的计算
Dist1 = (dist2(as.matrix(Data1), as.matrix(Data1)))^(1/2)
Dist2 = (dist2(as.matrix(Data2), as.matrix(Data2)))^(1/2)
代码语言:javascript复制##数据领域矩阵的构建
## First, set all the parameters:
K = 20; ##number of neighbors, must be greater than 1. usually (10~30)
alpha = 0.5; ##hyperparameter, usually (0.3~0.8)
T = 20; ###Number of Iterations, usually (10~50)
C = 2 ##number of cluster set by author
W1 = affinityMatrix(Dist1, K, alpha)
W2 = affinityMatrix(Dist2, K, alpha)
代码语言:javascript复制##载入标签
truelabel = c(matrix(1,100,1),matrix(2,100,1))
代码语言:javascript复制##可视化基于标签的聚类热图
displayClusters(W1, truelabel)
代码语言:javascript复制displayClusters(W2, truelabel)
代码语言:javascript复制##构建相似融合网络的数据
W = SNF(list(W1,W2), K, T)
代码语言:javascript复制##谱聚类的实现
labels = spectralClustering(W, C)
代码语言:javascript复制##计算互信息熵NMI。此处NMI范围[0,1],越大说明两者直接相关性越高,信息越一致。
NMI = calNMI(labels,truelabel)
代码语言:javascript复制##可视化谱聚类结果
###仅展示谱聚类结果标签
displayClustersWithHeatmap(W, labels)
代码语言:javascript复制###将其它分组信息加入展示
M_label=cbind(labels,truelabel)
colnames(M_label)=c("spectralClustering","TrueLabel")
M_label_colors=cbind("spectralClustering"=getColorsForGroups(M_label[,"spectralClustering"],
colors=c("blue","green")),"TrueLabel"=getColorsForGroups(M_label[,"TrueLabel"],
colors=c("orange","cyan")))
displayClustersWithHeatmap(W, labels, M_label_colors)
当然,我们还可以通过此包自定义的方法进行聚类数量的评估然后根据评估结果选择合适的簇数。此工具提供了两种评估算法,一种是启发式的特征值差值搜索(eigengap heuristic),其定义为若前k个特征值很小,且第k 1个特征值与前一个特征值相差比较大,则簇类个数选择k。这里的特征指的标准化的拉普拉斯矩阵的特征值,通过寻找前K个最小的特征并于K 1个存在很大差距。从而确定簇数为K;另一种是Rotation cost(未发现其相关的介绍)。下面直接看下代码:
代码语言:javascript复制estimateNumberOfClustersGivenGraph(W, NUMC=2:5)