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

节点可预测性——由精度矩阵直接得到每个症状的 R²

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

对标准化变量,把症状 i 对其余全部症状做线性回归,残差方差为 1/K_ii, 因此可预测性 R²_i = 1 − 1/K_ii(Haslbeck & Waldorp, 2018 的连续变量版本,此处为样本内值)。 代码同时用逐个回归验证这一等式。 输出解读:高可预测性的症状(如疲乏、情绪低落)被其邻居解释得多;可预测性低的症状 主要由网络外的因素决定,从网络内部”干预邻居”去改变它,预期效果有限。

Output输出图

节点可预测性——由精度矩阵直接得到每个症状的 R²:输出图 1

Output文字输出

情绪低落     R² = 0.417(逐个回归:0.417)
疲乏       R² = 0.401(逐个回归:0.401)
兴趣减退     R² = 0.316(逐个回归:0.316)
自我评价低    R² = 0.235(逐个回归:0.235)
注意力差     R² = 0.213(逐个回归:0.213)
睡眠问题     R² = 0.207(逐个回归:0.207)
自杀意念     R² = 0.189(逐个回归:0.189)
精神运动改变   R² = 0.068(逐个回归:0.068)
食欲改变     R² = 0.060(逐个回归:0.060)

Code代码

# 对标准化变量,把症状 i 对其余全部症状做线性回归,残差方差为 1/K_ii,
# 因此可预测性 R²_i = 1 − 1/K_ii(Haslbeck & Waldorp, 2018 的连续变量版本,此处为样本内值)。
# 代码同时用逐个回归验证这一等式。
# 输出解读:高可预测性的症状(如疲乏、情绪低落)被其邻居解释得多;可预测性低的症状
# 主要由网络外的因素决定,从网络内部"干预邻居"去改变它,预期效果有限。
import numpy as np

TRUE_EDGES = {(0, 1): .40, (1, 5): .30, (1, 8): .25, (5, 8): .25, (2, 3): .40,
              (3, 6): .30, (3, 7): .20, (3, 4): .20, (0, 3): .20, (6, 7): .20}
names = ["兴趣减退", "情绪低落", "睡眠问题", "疲乏", "食欲改变",
         "自我评价低", "注意力差", "精神运动改变", "自杀意念"]
K0 = np.eye(9)
for (i, j), w in TRUE_EDGES.items():
    K0[i, j] = K0[j, i] = -w
rng = np.random.default_rng(0)
X = rng.multivariate_normal(np.zeros(9), np.linalg.inv(K0), size=500)   # 模拟数据,仅用于演示方法
Z = (X - X.mean(0)) / X.std(0)

R = np.corrcoef(Z, rowvar=False)
K = np.linalg.inv(R)
r2_formula = 1 - 1 / np.diag(K)
r2_reg = []
for i in range(9):
    A = np.delete(Z, i, axis=1)
    b, *_ = np.linalg.lstsq(A, Z[:, i], rcond=None)
    r2_reg.append(1 - np.var(Z[:, i] - A @ b) / np.var(Z[:, i]))
for i in np.argsort(r2_formula)[::-1]:
    print(f"{names[i]:<8} R² = {r2_formula[i]:.3f}(逐个回归:{r2_reg[i]:.3f})")

import networkx as nx
import matplotlib.pyplot as plt
abbr = ["Anh", "Sad", "Slp", "Fat", "App", "Gui", "Con", "Mot", "Sui"]
dk = np.sqrt(np.diag(K)); Pc = -K / np.outer(dk, dk); np.fill_diagonal(Pc, 0)
fig, axes = plt.subplots(1, 2, figsize=(10, 4.2))
pos = nx.circular_layout(range(9))
G = nx.Graph(); G.add_nodes_from(range(9))
G.add_weighted_edges_from((i, j, Pc[i, j]) for i in range(9) for j in range(i + 1, 9) if abs(Pc[i, j]) > 0.1)
ws = [d["weight"] for _, _, d in G.edges(data=True)]
nx.draw_networkx_edges(G, pos, ax=axes[0], width=[max(1.0, 8 * abs(w)) for w in ws],
                       edge_color=["#4477AA" if w > 0 else "#CC6677" for w in ws])
nodes = nx.draw_networkx_nodes(G, pos, ax=axes[0], node_size=200, node_color=r2_formula, cmap="Oranges",
                               vmin=0, vmax=r2_formula.max(), edgecolors="k")
nx.draw_networkx_labels(G, {k: 1.22 * v for k, v in pos.items()}, labels=dict(enumerate(abbr)), ax=axes[0], font_size=9)
fig.colorbar(nodes, ax=axes[0], shrink=0.7, label="$R^2$")
axes[0].set_title("Node colour = predictability (|pcor| > 0.1 shown)", fontsize=10)
axes[0].set_aspect("equal"); axes[0].axis("off"); axes[0].set_xlim(-1.45, 1.45); axes[0].set_ylim(-1.45, 1.45)
order = np.argsort(r2_formula)
axes[1].barh(np.arange(9), r2_formula[order], color="#4477AA", label="$1 - 1/K_{ii}$")
axes[1].scatter(np.array(r2_reg)[order], np.arange(9), marker="|", color="k", s=200, zorder=3, label="per-node regression")
axes[1].set_yticks(range(9), [abbr[i] for i in order])
axes[1].set(xlabel="Node predictability $R^2$", title="Precision matrix gives $R^2$ directly")
axes[1].legend(fontsize=8, loc="lower right")
plt.tight_layout(); plt.show()

Run在浏览器中运行

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

节点可预测性——由精度矩阵直接得到每个症状的 R²