06 SNF 相似性网络与亚型识别
SNF（Similarity Network Fusion）的思路和 MOFA 完全不同：它不做矩阵分解，而是为每一层组学数据分别构建样本相似性网络，然后把这些网络融合成一个统一的网络。融合后的网络可以用谱聚类（spectral clustering）来识别亚型。
算法直觉
- 对每一层组学数据，计算样本间的距离矩阵（通常用欧氏距离）
- 把距离矩阵转换成相似性矩阵（用 scaled exponential kernel）
- 构建 K 近邻图（KNN graph）
- 迭代融合：每一层的网络向其他层"扩散"信息，最终收敛到一个融合网络
融合的核心思想是：如果两个样本在多层数据中都相似，它们在融合网络中的连接会被加强；如果只在一层相似，连接会被削弱。
基本流程
library(SNFtool)
# 输入：每层数据矩阵，行是样本，列是特征
rna_t <- t(rna_final)
meth_t <- t(meth_final)
prot_t <- t(prot_final)
# 参数
K <- 20 # 近邻数
alpha <- 0.5 # 超参数
iter <- 20 # 迭代次数
# 1. 计算距离矩阵
dist_rna <- dist2(as.matrix(rna_t), as.matrix(rna_t))
dist_meth <- dist2(as.matrix(meth_t), as.matrix(meth_t))
dist_prot <- dist2(as.matrix(prot_t), as.matrix(prot_t))
# 2. 构建相似性网络
W_rna <- affinityMatrix(dist_rna, K, alpha)
W_meth <- affinityMatrix(dist_meth, K, alpha)
W_prot <- affinityMatrix(dist_prot, K, alpha)
# 3. 融合网络
W_fused <- SNF(list(W_rna, W_meth, W_prot), K, iter)
谱聚类确定亚型
# 用 eigen-gap 方法估计最优聚类数
C_best <- estimateNumberOfClustersGivenGraph(W_fused, NUMC = 2:8)
cat("推荐聚类数:", C_best$`Eigen-gap best`, "\n")
# 谱聚类
clusters <- spectralClustering(W_fused, K = C_best$`Eigen-gap best`)
table(clusters)
# 可视化融合网络
library(pheatmap)
pheatmap(W_fused,
annotation_row = data.frame(Cluster = factor(clusters)),
show_rownames = FALSE, show_colnames = FALSE,
main = "Fused similarity network")
与 Consensus Clustering 比较
SNF 和 consensus clustering（如 ConsensusClusterPlus）都能做亚型发现，但思路不同：
| 方面 | SNF | Consensus Clustering |
|---|---|---|
| 输入 | 多层数据分别构建网络 | 通常用拼接后的单矩阵 |
| 整合方式 | 网络融合（中期整合） | 重采样 + 聚类稳定性（早期整合） |
| 对噪声层的鲁棒性 | 较好，噪声层贡献会被削弱 | 较差，噪声层会污染拼接矩阵 |
| 缺失值 | 不支持 | 不支持 |
| 聚类数选择 | eigen-gap | CDF / delta area |
实际项目中可以两种都跑，看结果是否一致。如果一致，结论更可信。
library(ConsensusClusterPlus)
# 拼接矩阵
combined <- rbind(rna_final, meth_final, prot_final)
# Consensus clustering
cc_results <- ConsensusClusterPlus(combined,
maxK = 8,
reps = 500,
pItem = 0.8,
pFeature = 1,
clusterAlg = "hc",
distance = "pearson",
seed = 42,
plot = "pdf")
亚型与临床表型关联
# 把聚类结果和临床信息合并
cluster_df <- data.frame(
sample = rownames(rna_t),
cluster = factor(clusters),
subtype = col_data$subtype,
stage = col_data$stage
)
# Fisher 检验：聚类 vs 已知亚型
fisher.test(table(cluster_df$cluster, cluster_df$subtype))
# 生存分析
library(survival)
library(survminer)
surv_obj <- Surv(col_data$os_time, col_data$os_event)
fit <- survfit(surv_obj ~ cluster_df$cluster)
ggsurvplot(fit, data = cluster_df, pval = TRUE,
palette = "jco", risk.table = TRUE)
参数敏感性
SNF 对 K（近邻数）比较敏感。建议在 K = 10-30 范围内扫描，看聚类结果是否稳定。alpha 通常固定为 0.5，迭代次数 20 次足够收敛。
常见坑
坑 1：K 值敏感但没做参数扫描
SNF 对近邻数 K 非常敏感。K = 10 和 K = 30 可能给出完全不同的聚类。很多人只跑一次 K = 20 就报告结果。正确做法是在 K = 10-30 范围内扫描，看聚类标签的 NMI（normalized mutual information）是否稳定。如果不同 K 结果差异大，结论不可信。
坑 2：样本量太小时网络融合不可靠
SNF 构建 KNN 图时，如果样本数 < 30 而 K = 20，几乎每个样本都是彼此的邻居，网络失去区分力。小样本场景建议把 K 降到 sample_size / 5 以下，或者换用 MOFA2（对小样本更友好）。
坑 3：某一层数据质量差污染融合结果
SNF 会把噪声层的信息也融合进去。如果蛋白组只有 500 个特征且缺失率高，它的相似性矩阵很嘈杂，融合后会拖累其他层的信号。建议先单独检查每层的 silhouette score，质量太差的层考虑排除。
坑 4：用 eigen-gap 选聚类数时只看一个指标
estimateNumberOfClustersGivenGraph 返回多个指标（eigen-gap、rotation cost 等），它们可能给出不同的最优 K。不要只看一个。如果多个指标一致推荐 K = 3，结果更可信；如果分歧大，说明数据中没有清晰的亚型结构。
下一步
接着深入：
- 07 机器学习预测模型 — 从亚型发现转向临床预测
- 08 整合结果的生物学解读 — 把聚类标签翻译成生物学亚型
横向延伸：
- 05 MOFA2 因子分析实战 — 另一种整合方法，适合小样本和部分缺失
- bulk RNA-seq 07 聚类与亚型 — 单组学聚类方法的对比
参考资源
离线资料下载
手册 HTML / PDF 已在后台预 生成，点击后直接下载网站静态资源。