生物医学网络科学Biomedical Network Science
Example07 / 10

网络比较检验（NCT）的置换逻辑

所属专题
症状网络与网络精神病学
下载
Notebook · .py

两组被试（模拟数据，仅用于演示方法）：B 组在”睡眠问题–疲乏”上的真实偏相关更强，其余相同。 统计量：M = 两组边权差的最大绝对值（网络结构不变性），S = 两组全局强度之差的绝对值（全局强度不变性）。 零分布：把两组被试合并后随机重分组（保持组大小），重复 1000 次重新计算 M 与 S。 输出解读：M 的置换 p 值小，说明至少一条边不同；再逐边检验并做多重比较校正才能定位是哪条边。 S 的 p 值不显著：一条边变强并不一定使全局强度显著改变，三类检验回答的是不同问题。 简化：这里用未正则化偏相关以便快速演示；NetworkComparisonTest 默认在每次置换中重新运行 EBICglasso。

Output输出图

网络比较检验（NCT）的置换逻辑：输出图 1

Output文字输出

网络结构不变性：M = 0.285，置换 p = 0.001
全局强度不变性：S = 0.124，置换 p = 0.491
Holm 校正后显著不同的边：2-3，|Δw| = 0.285，p_holm = 0.015

Code代码

# 两组被试（模拟数据，仅用于演示方法）：B 组在"睡眠问题–疲乏"上的真实偏相关更强，其余相同。
# 统计量：M = 两组边权差的最大绝对值（网络结构不变性），S = 两组全局强度之差的绝对值（全局强度不变性）。
# 零分布：把两组被试合并后随机重分组（保持组大小），重复 1000 次重新计算 M 与 S。
# 输出解读：M 的置换 p 值小，说明至少一条边不同；再逐边检验并做多重比较校正才能定位是哪条边。
# S 的 p 值不显著：一条边变强并不一定使全局强度显著改变，三类检验回答的是不同问题。
# 简化：这里用未正则化偏相关以便快速演示；NetworkComparisonTest 默认在每次置换中重新运行 EBICglasso。
import numpy as np

def make_K(extra=0.0):
    K = np.eye(6)
    for (i, j), w in {(0, 1): .3, (1, 2): .25, (2, 3): .3, (3, 4): .2, (4, 5): .3, (0, 5): .2}.items():
        K[i, j] = K[j, i] = -w
    K[2, 3] = K[3, 2] = -(0.3 + extra)   # 节点 2 = 睡眠问题，节点 3 = 疲乏
    return K

def pcor(X):
    K = np.linalg.inv(np.cov(X, rowvar=False))
    d = np.sqrt(np.diag(K))
    P = -K / np.outer(d, d)
    np.fill_diagonal(P, 0)
    return P

rng = np.random.default_rng(7)
nA, nB = 400, 400
XA = rng.multivariate_normal(np.zeros(6), np.linalg.inv(make_K()), size=nA)
XB = rng.multivariate_normal(np.zeros(6), np.linalg.inv(make_K(0.25)), size=nB)
iu = np.triu_indices(6, 1)

def stats(XA, XB):
    D = pcor(XA)[iu] - pcor(XB)[iu]
    S = abs(np.abs(pcor(XA)[iu]).sum() - np.abs(pcor(XB)[iu]).sum())
    return np.max(np.abs(D)), S, np.abs(D)

M0, S0, D0 = stats(XA, XB)
Xall = np.vstack([XA, XB])
perm_M, perm_S, perm_D = [], [], []
for _ in range(1000):
    idx = rng.permutation(nA + nB)
    m, s, d = stats(Xall[idx[:nA]], Xall[idx[nA:]])
    perm_M.append(m); perm_S.append(s); perm_D.append(d)
perm_D = np.array(perm_D)
p_M = (1 + np.sum(np.array(perm_M) >= M0)) / 1001
p_S = (1 + np.sum(np.array(perm_S) >= S0)) / 1001
print(f"网络结构不变性：M = {M0:.3f}，置换 p = {p_M:.3f}")
print(f"全局强度不变性：S = {S0:.3f}，置换 p = {p_S:.3f}")
p_edge = (1 + (perm_D >= D0).sum(axis=0)) / 1001
m = len(p_edge)                          # Holm 校正
order = np.argsort(p_edge)
p_holm = np.empty(m)
p_holm[order] = np.minimum(1, np.maximum.accumulate((m - np.arange(m)) * p_edge[order]))
for k, (i, j) in enumerate(zip(*iu)):
    if p_holm[k] < 0.05:
        print(f"Holm 校正后显著不同的边：{i}-{j}，|Δw| = {D0[k]:.3f}，p_holm = {p_holm[k]:.3f}")

import matplotlib.pyplot as plt
fig, axes = plt.subplots(1, 2, figsize=(10, 4))
axes[0].hist(perm_M, bins=30, color="#BBBBBB", edgecolor="k", lw=0.4, label="permutations")
axes[0].axvline(M0, color="crimson", lw=2, label=f"observed M (p = {p_M:.3f})")
axes[0].set(xlabel="Max |edge difference| M", ylabel="Count", title="Network structure invariance test")
axes[0].legend(fontsize=8)
x = np.arange(len(D0))
axes[1].bar(x, D0, color=np.where(p_holm < 0.05, "crimson", "#4477AA"))
axes[1].scatter(x, np.quantile(perm_D, 0.95, axis=0), marker="_", s=300, color="k", label="permutation 95th pct")
axes[1].set_xticks(x, [f"{i}-{j}" for i, j in zip(*iu)], fontsize=8)
axes[1].set(xlabel="Edge (node 2 = sleep, 3 = fatigue)", ylabel="|pcor_A - pcor_B|", title="Edge-wise tests (red: Holm p < 0.05)")
axes[1].legend(fontsize=8)
plt.tight_layout(); plt.show()

Run在浏览器中运行

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

网络比较检验（NCT）的置换逻辑