完整案例：从表达矩阵到可检验的结论

本节用一个模拟的共表达数据集，把前面各节串成一次完整分析。数据由代码生成，含已知答案（植入的模块与疾病基因），因此可以检查每一步是否把答案找回来了；真实分析的步骤相同，只是没有答案可对。完整代码见示例“完整案例”。

研究问题：一组候选疾病基因在共表达网络中是否比度相近的随机基因更聚集？若是，它们集中在哪个模块，模块的枢纽基因是哪些？

数据：200 个基因 × 60 个样本。160 个基因分属 4 个各含 40 个基因的模块，每个模块由一个潜在因子驱动（因子载荷 0.5–0.9）；另有 40 个噪声基因，其中 5 个几乎不随样本变化。各基因的表达尺度不同。疾病基因集含 12 个基因：8 个来自模块 0，4 个随机取自其他基因。

第 1 步，质控：去除方差极低的基因（表达几乎不变的基因，其相关系数主要由噪声决定）。去除 5 个，保留 195 个。

第 2 步，估计边并控制假阳性：计算 Spearman 相关（对异常值和非线性单调关系稳健，且不受表达尺度影响）。用 Fisher $z$ 变换求每个相关系数的 $p$ 值——对 Spearman 相关，$\operatorname{arctanh}(r_s)$ 的方差约为 $1.06/(n-3)$（Fieller 等 1957）——再用 Benjamini–Hochberg 方法在全部约 1.9 万个基因对上把错误发现率（FDR）控制在 0.01，只保留正相关。得到 2031 条边，密度约 0.11；34 个基因没有任何边，其中 33 个是噪声基因。这一步把“阈值取多少”转化为“容忍多大比例的假边”，阈值的含义更清楚（另见阈值、加权与稀疏化）。

第 3 步，模块识别：在去除孤立点后的加权网络上运行 Louvain（固定种子），得到 4 个模块，$Q\approx0.75$；与植入模块对照，调整兰德指数（ARI）约 0.98。真实数据没有“植入模块”，此时应当报告多次运行的稳定性，并用独立数据或功能注释检验模块（见社团结果的验证与解读）。

第 4 步，枢纽基因：在每个模块内部计算加权度（模块内连接度，与 WGCNA 中 intramodular connectivity 的思想一致，见 WGCNA），列出每个模块中最高者作为候选枢纽。

第 5 步，零模型检验：疾病基因集中有 11 个在网络中（1 个随机选入的噪声基因成为孤立点而被去除）。用两种统计量衡量“聚集”：诱导子图最大连通分支（LCC）的大小，以及集合内部的边数。零模型是度匹配的随机基因集：把网络节点按度分为 10 个分位区间，每个疾病基因替换为同区间内的随机基因，重复 1000 次。结果：LCC 观测值 9，随机期望 4.5 ± 1.0，$z\approx4.7$；内部边数观测值 24，随机期望 10.1 ± 2.7，$z\approx5.1$；两者的经验 $p$ 值均约为 0.003。

第 6 步，敏感性分析：把 FDR 换成 0.001、0.05、0.1 重建网络并重复检验。内部边数在三种设置下都显著（$p$ 约 0.002–0.018）；LCC 在 FDR 0.001 与 0.05 下显著，在 FDR 0.1 时不再显著（观测 9，随机期望 7.0，$p\approx0.21$）。原因是网络变稠密后，随机基因集也很容易连成大分支，LCC 这个统计量趋于饱和。这提示：统计量的选择与网络密度会相互作用，应当预先指定主要统计量，并报告全部敏感性结果，而不是挑选显著的那个。

结论的写法：“在以 Spearman 相关、FDR 0.01 构建的共表达网络中，候选疾病基因的内部连接显著多于度匹配的随机基因集（$z\approx5.1$，经验 $p\approx0.003$），且对 FDR 水平稳健；以 LCC 衡量时，结论在较稠密的网络中不再成立。”这句话限定了网络构建方式、零模型和统计量，比“这些基因构成致病模块”更准确，也更容易被他人检验。