生物医学网络科学Biomedical Network Science
Example13 / 13

相似性网络融合(SNF)整合两类组学数据

所属专题
生物分子网络与拓扑特征
下载
Notebook · .py

模拟 90 名患者分属 3 个亚型;“表达”数据只能区分亚型 1 与 {2,3},“甲基化”数据只能区分亚型 3 与 {1,2}。 对每类数据建患者相似性矩阵(缩放指数核),再按 Wang 等(2014)的思路迭代: P_v ← S_v · (其他数据类型的 P 的平均) · S_vᵀ,其中 S_v 只保留每个患者的 K 个近邻。 最后对融合矩阵做谱聚类(拉普拉斯特征向量 + k-means),与真实亚型比较调整兰德指数(ARI)。 输出解读:单一数据类型只能分出 2 组,融合后可分出 3 个亚型(ARI 明显更高)。

Output输出图

相似性网络融合(SNF)整合两类组学数据:输出图 1

Output文字输出

只用表达:ARI = 0.47
只用甲基化:ARI = 0.31
SNF 融合:ARI = 0.83

Code代码

# 模拟 90 名患者分属 3 个亚型;“表达”数据只能区分亚型 1 与 {2,3},“甲基化”数据只能区分亚型 3 与 {1,2}。
# 对每类数据建患者相似性矩阵(缩放指数核),再按 Wang 等(2014)的思路迭代:
#   P_v ← S_v · (其他数据类型的 P 的平均) · S_vᵀ,其中 S_v 只保留每个患者的 K 个近邻。
# 最后对融合矩阵做谱聚类(拉普拉斯特征向量 + k-means),与真实亚型比较调整兰德指数(ARI)。
# 输出解读:单一数据类型只能分出 2 组,融合后可分出 3 个亚型(ARI 明显更高)。
import numpy as np
from scipy.cluster.vq import kmeans2
from scipy.spatial.distance import cdist
from scipy.special import comb

rng = np.random.default_rng(4)
labels = np.repeat([0, 1, 2], 30)                   # 模拟数据,仅用于演示方法
expr = rng.normal(size=(90, 50)); expr[labels == 0, :10] += 1.5
meth = rng.normal(size=(90, 40)); meth[labels == 2, :10] += 1.5

def affinity(X, K=15, mu=0.5):
    D = cdist(X, X); Ds = np.sort(D, 1)[:, 1:K + 1].mean(1)
    eps = (Ds[:, None] + Ds[None, :] + D) / 3
    return np.exp(-D**2 / (mu * eps**2))

def full_kernel(W):
    P = W / (2 * (W.sum(1, keepdims=True) - np.diag(W)[:, None])); np.fill_diagonal(P, 0.5); return P

def knn_kernel(W, K=15):
    S = np.zeros_like(W)
    for i in range(len(W)):
        nb = np.argsort(-W[i])[:K]; S[i, nb] = W[i, nb] / W[i, nb].sum()
    return S

def snf(Ws, K=15, t=20):
    Ps = [full_kernel(W) for W in Ws]; Ss = [knn_kernel(W, K) for W in Ws]
    for _ in range(t):
        Ps = [Ss[v] @ np.mean([Ps[u] for u in range(len(Ws)) if u != v], 0) @ Ss[v].T for v in range(len(Ws))]
        Ps = [(P + P.T) / 2 for P in Ps]
    return np.mean(Ps, 0)

def spectral(W, k, seed=0):
    d = W.sum(1); L = np.eye(len(W)) - W / np.sqrt(np.outer(d, d))
    _, vec = np.linalg.eigh(L); U = vec[:, :k]; U /= np.linalg.norm(U, axis=1, keepdims=True)
    return kmeans2(U, k, seed=seed, minit="++")[1]

def ari(a, b):
    ct = np.array([[np.sum((a == i) & (b == j)) for j in np.unique(b)] for i in np.unique(a)])
    s = comb(ct, 2).sum(); sa = comb(ct.sum(1), 2).sum(); sb = comb(ct.sum(0), 2).sum(); tot = comb(len(a), 2)
    exp = sa * sb / tot
    return (s - exp) / (0.5 * (sa + sb) - exp)

W_e, W_m = affinity(expr), affinity(meth)
W_f = snf([W_e, W_m])
for name, W in [("只用表达", W_e), ("只用甲基化", W_m), ("SNF 融合", W_f)]:
    print(f"{name}:ARI = {ari(labels, spectral(W, 3)):.2f}")

import matplotlib.pyplot as plt
fig, axes = plt.subplots(1, 3, figsize=(11, 3.8))
for ax, W, title in zip(axes, [W_e, W_m, W_f], ["expression only", "methylation only", "SNF fused"]):
    M = W.copy(); np.fill_diagonal(M, np.nan)              # 患者按真实亚型排序,去掉对角线以便配色
    ax.imshow(M, cmap="viridis", vmax=np.nanpercentile(M, 98), interpolation="nearest")
    for b in (29.5, 59.5):
        ax.axhline(b, color="w", lw=0.8); ax.axvline(b, color="w", lw=0.8)
    ax.set_title(title); ax.set_xticks([15, 45, 75]); ax.set_xticklabels(["S1", "S2", "S3"])
    ax.set_yticks([15, 45, 75]); ax.set_yticklabels(["S1", "S2", "S3"])
fig.suptitle("Patient similarity ordered by true subtype: fusion reveals 3 blocks")
plt.tight_layout(); plt.show()

Run在浏览器中运行

首次加载 Python 环境约需 20–40 秒;可以修改参数后重新运行。

相似性网络融合(SNF)整合两类组学数据