简介：面向神经科学领域研究者与脑影像分析初学者，聚焦结构协变网络（SCN）的构建与分析，解决从MRI/DTI等原始影像到网络指标计算的全流程问题。资源覆盖数据获取、图像预处理、区域分割、特征提取、协方差计算、网络构建及网络分析等关键环节，可用于模块检测、中心性测量、小世界性分析等场景，在认知功能研究、阿尔茨海默病与自闭症等疾病诊断、脑损伤恢复评估中均有应用价值。压缩包共239个文件，以MATLAB脚本（.m）为主，辅以矩阵数据（.mat）、医学图像模板（.nii/.hdr/.img）及示例图谱（.bmp/.jpg），整体约19.19MB，既有标准模板文件也有可直接修改运行的脚本，便于用户按需组合或二次开发。已有1109人学习下载，包内目录结构清晰，适合有一定MATLAB基础、希望系统掌握SCN构建流程的研究者快速上手并通过示例数据复现各步骤结果。

1. 结构协变网络：别把个体网络的思路套到群体数据上

拿到一批 T1 结构像，想在脑区之间建一张“谁和谁一起变大变小”的网络，这个需求对应的方法就是结构协变网络（Structural Covariance Network, SCN）。NeuroimageTools-main 背后做的事，是把 T1 分割结果整理成一张“被试 × 脑区”的形态指标矩阵，再在组水平上算出脑区之间的协变关系——它回答的不是“某个人脑内哪里连着哪里”，而是“这一组人里，哪些区域的灰质体积或皮层厚度会同步升高、同步下降”。

这个区别是很多人翻车的起点。功能连接（fMRI）可以给每个被试单独出一张网，结构协变网络却做不到：一个被试只有一组数值，算不出协变。必须带着一整组被试，比如 60 个患者、60 个对照，才能得到一张有统计意义的组水平网络。适合做这件事的人，是手里有 T1 和分割结果、想比较两组人群网络拓扑差异的研究者；不适合的人，是期待“输一个被试就能出图”的临床工具使用者。下面我按自己跑通这条流程的顺序，把构建 SCN 的每一步、每个参数和每个坑写清楚。

2. 协变网络的两个数学前提：组水平共变与形态指标选型

2.1 组水平设计矩阵：SCN 的根本定义

结构协变网络的计算对象不是脑图像，而是一个矩阵。假设有 N 个被试，每个被试的 T1 经过分割后，可以得到 M 个脑区的形态学指标（灰质体积、皮层厚度、表面积等）。把这些数据排列成设计矩阵 X，形状是 N×M，每一列是一个脑区在所有被试上的指标取值。

SCN 的构建，就是对 X 的每一对列做相关性分析，得到一个 M×M 的协变矩阵。以灰质体积为例：如果左海马和左杏仁核两列数值在被试之间同步升降，Pearson 相关系数就高，对应边权重就大。整个过程没有用到任何个体内部的连接信息——没有纤维束，没有 BOLD 时间序列，没有功能激活。它描述的是群体层面的形态学共变模式，这也是“协变”二字的数学含义。

这里要和另一种做法划清界限：用单个被试 T1 的体素形态相似性去构建个体结构网络（有些人叫它 T1w 个体形态网络）。这种思路看起来更诱人，但你的输入是一组被试而不是一个人，且它对平滑核、体素分辨率、配准质量极度敏感，很多实验室重复出来的结果不稳定。我一般不建议新人在这个方向上耗费时间：SCN 的组水平逻辑更稳健、统计检验更成熟、审稿人接受度也更高。如果你明确要“组间比较”，直接走组水平 SCN 是风险最低的路径。

2.2 形态学指标选哪个：GMV、皮层厚度、表面积的取舍

指标选择直接决定网络生物学意义，这一步没有“默认最优”，只有“最匹配你的问题”。下表是我在实际项目中的选择逻辑：

经验是：回答“哪里萎缩了”选皮层厚度；回答“哪个系统退化重”选 GMV；想证明结果不是指标偶然，就用两种指标各跑一遍看拓扑一致性。皮层下结构通常用体积指标（来自分割标签的统计表），皮层区域用厚度或体积都可以。

2.3 模板与节点粒度：68 节点还是 148 节点，被试数够不够

节点数由分区模板决定。我常用的是 DK 模板（68 个皮层区），它节点数适中，空间平滑后信噪比好，跨研究可比性高。另一套 148 节点的模板更精细，能捕捉更细的模块边界，但对分割误差更敏感，且对样本量要求更高。

这里有一条硬约束：如果你打算用偏相关构建网络（控制其他脑区影响），矩阵求逆要求被试数显著大于节点数。我一般要求 N/M 至少大于 3，最好达到 5。也就是说：68 节点至少 200 被试起步，148 节点至少 400 被试。达不到时，要么改用 Pearson 相关，要么用收缩估计（Ledoit-Wolf），不能硬算偏相关，否则矩阵非正定，后面所有图论指标都会变成玄学。

节点粒度的另一个影响是模块检测：节点越密，模块划分越碎，跨研究重复性越差。如果你后续想做模块化分析，68 节点比 148 节点更稳；想做 hub 节点定位，148 节点能给出更精确的坐标。

3. 构建流程落地：从 T1 分割到协变矩阵的全过程

3.1 数据准备与分割质量控制

SCN 的起点不是网络，而是分割质量。这一步决定后面所有结果的可靠性，比任何参数调优都重要。

常见做法是在 T1 分割流水线跑完后，先做两项检查：一是自动质量指标，重点看分割包生成的脑掩膜与原始 T1 的配准情况，以及表面重建的拓扑缺陷数量；二是人工目检，把每个被试的脑掩膜叠在 T1 上，逐层排查明显的脑组织缺失或颅骨残留。我曾见过一个数据集，某被试分割时把双侧颞叶外缘切掉了一部分，如果只看自动指标根本发现不了，但它让那个被试的颞叶 GMV 直接少了 15%，进入矩阵后形成一条虚假的强协变边。

把数据从分割结果提取成表格这一步，可以用分割工具自带的统计脚本，把每个 ROI 的体积和厚度汇总成一张表，然后用下面这段代码做初步筛查：

# 对每个被试，计算皮层各 ROI 的平均厚度（示意流程） # 输入: 某被试的皮层表面厚度图与分区标签 # 输出: roi_thickness.csv，行列分别为 ROI 与被试ID for sub in subject_list; do mri_segstats --subject "$sub" \ --i "$sub/surf/lh.thickness" \ --annot "$sub/label/lh.DK.aparc.annot" \ --sum "$sub/stats/lh.DK.thickness.stats" # 右半球同理，这里省略 done

这段命令做的事是：对左半球表面，把厚度值按 DK 模板的每个 ROI 取均值，输出成统计文件。参数上，--i指定输入图（厚度图），--annot指定表面分区标签，--sum指定输出路径。实际环境里命令名会随工具版本不同，但“按标签统计表面值”这个思路是通用的。

跑完后把左右半球和皮层下体积合并成宽表：行是被试，列是 ROI 名称。合并时最容易出的错是标签错位，我放在第 5 章单独讲。质量控制阶段还要看每个 ROI 的分布：超过均值 3 个标准差的被试标记出来，回到原始 T1 重新目检，不轻易删除，但一定要记录。

3.2 清洗、协变量回归与设计矩阵构造

拿到原始指标表后，第一步不是算相关，而是去掉你不想让网络编码的变异。结构协变网络最怕两件事：年龄效应和颅内总体积（ICV）效应。年龄与几乎所有脑区的灰质体积都强相关，ICV 与体积指标同样强相关——不回归这两个变量，网络会被“全局衰老”或“脑袋大小”主导，得到一张全脑几乎全连通的图。

我的做法是：先做线性回归，把年龄、性别、ICV 以及扫描站点（如果多中心）从每个 ROI 的指标中剔除，用残差作为后续计算输入。这一步用 Python 可以很干净地完成：

import pandas as pd import numpy as np from sklearn.linear_model import LinearRegression # 读取宽表：行为被试，列为 ROI，另有协变量列 df = pd.read_csv("roi_wide.csv") # 包含 age, sex, icv 与 68 个 ROI roi_cols = [c for c in df.columns if c.startswith("roi_")] # 构造协变量矩阵: 年龄 + 性别 + 颅内总体积 # 性别已编码为 0/1; ICV 单位建议用 mm^3 并做 log 转换后使用 X_cov = df[["age", "sex", "icv_log"]].values resid = pd.DataFrame(index=df.index, columns=roi_cols) for roi in roi_cols: reg = LinearRegression().fit(X_cov, df[roi].values) resid[roi] = df[roi].values - reg.predict(X_cov) resid.to_csv("roi_resid.csv", index=False)

关键点：回归是逐 ROI 做的，每个 ROI 拥有独立的回归系数，因为年龄对海马和初级视觉皮层的效应强度完全不同。icv_log是对 ICV 取对数，减少极端值影响。回归后的残差应该近似零均值、无趋势，可以画一两个 ROI 的残差分布确认正态性。

如果你用的是皮层厚度，ICV 是否回归存在争议——厚度与 ICV 相关性本来就弱，回归与否影响不大；但 GMV 必须回归。我的习惯是：GMV 一律回归 ICV；厚度在结果里同时报告“回归 ICV”和“不回归 ICV”两种版本，证明结论不依赖这个决定。

3.3 相关矩阵与偏相关矩阵：从 Pearson 到 Ledoit-Wolf 收缩

设计矩阵干净之后，就到了构建网络矩阵的核心步骤。这里有三条技术路线，难度和稳健性递进：

第一，Pearson 相关：最简单直接，每个 ROI 对算相关系数，得到 M×M 矩阵。缺点是没有控制其他脑区影响，网络中的边可能只是间接相关。

第二，偏相关：在控制其他所有 ROI 的条件下，计算两两相关。这更接近“直接连接”的语义，但要求协方差矩阵可逆，样本量不足时结果极不稳定。

第三，收缩偏相关：用 Ledoit-Wolf 方法对协方差矩阵做收缩估计，在不牺牲太多精度的前提下保证矩阵正定。这是我在样本量不太充足时的默认选择：

import numpy as np import pandas as pd from sklearn.covariance import LedoitWolf resid = pd.read_csv("roi_resid.csv") X = resid.values # shape: (n_subjects, n_rois) # Ledoit-Wolf 收缩协方差估计 lw = LedoitWolf().fit(X) cov = lw.covariance_ # 由精度矩阵(逆协方差)计算偏相关 prec = np.linalg.inv(cov) d = np.sqrt(np.diag(prec)) partial_corr = prec / np.outer(d, d) np.fill_diagonal(partial_corr, 0) # 保存网络矩阵: 行/列均为 ROI pd.DataFrame(partial_corr, index=resid.columns, columns=resid.columns).to_csv("scn_partial_corr.csv")

LedoitWolf()会自动估计收缩强度——当样本量不足时它自动加强收缩，样本充足时收缩量趋近于零，可以理解为自动化的“保守模式”。np.outer(d, d)是在做矩阵归一化，把精度矩阵转换成偏相关矩阵，这一步不能省。

参数层面的建议：如果 N/M 大于 5，直接用普通偏相关即可；如果 N/M 在 3 到 5，用 Ledoit-Wolf；小于 3，建议回归到 Pearson 相关，或者考虑减少节点数。网络矩阵构建完成后，下一步是阈值化，这直接决定网络拓扑指标的长相。

4. 阈值化与拓扑参数：密度阈值、负连接与加权网络

4.1 为什么不要用 r>0.3 做绝对阈值

很多初学者拿到相关矩阵的第一反应是“把相关系数大于 0.3 的边留下，其余删掉”。这个做法在结构协变网络里基本不可取。

原因有二。第一，相关系数的绝对值大小受样本量影响：同一个真实相关强度，在大样本里估计值更稳定，在小样本里波动更大。两组被试样本量不同时，用绝对阈值会得到密度差异悬殊的网络，后续比较“聚类系数谁更高”时，你比较的其实是密度而不是拓扑。第二，SCN 的相关系数普遍不高，健康人群的脑区间协变强度通常在 0.2–0.6 之间，r>0.3 的阈值会丢掉大量真实的弱连接，让网络只剩最显眼的几条边。

我的经验是：网络比较的基准应该是“密度”，不是“阈值”。密度（或称稀疏度）定义为保留的边数占可能边数的比例，一般取 5%–40% 之间的几个水平，例如 5%、10%、15%、20%、25%、30%。每个密度水平下保留最强的那一批边，然后分别计算拓扑指标，最后看指标跨密度水平的一致性。这样两组被试的网络密度完全相同，比较才有意义。

4.2 密度阈值实现：多阈值稳健性检验

密度阈值化的代码比听起来简单，核心是对相关矩阵的上三角取绝对值排序，按密度截断：

import numpy as np def threshold_by_density(mat, density): """ 按密度阈值化相关/偏相关矩阵 mat: 对称矩阵, 对角线为0 density: 保留边比例, 如 0.1 表示保留最强者10% 返回: 加权邻接矩阵(不保留的边置0) """ mat = mat.copy() np.fill_diagonal(mat, 0) # 消除自环 tri = np.triu_indices_from(mat, k=1) # 取上三角避免重复 # 按绝对值排序, 取第 (density * n_possible_edges) 大的值作为阈值 vals = np.abs(mat[tri]) n_possible = len(vals) n_keep = int(np.ceil(n_possible * density)) kth = np.sort(vals)[::-1][n_keep - 1] # 低于阈值的边置0, 其余保留原始权重 mat[mat < kth] = 0 return mat # 对偏相关矩阵做 10% 密度阈值 for density in [0.05, 0.10, 0.15, 0.20, 0.25, 0.30]: net = threshold_by_density(partial_corr, density) # 在这里计算聚类系数、特征路径长度、全局效率等指标

这段代码的逻辑是：先取矩阵上三角的所有边的绝对值，排序后找到“第 n_keep 大”的值作为阈值，低于阈值的边全部置零，保留的边保留原始权重。参数density是你指定的网络稀疏度，越小网络越稀疏。注意最后一步用的是mat[mat < kth] = 0，因为mat是实值相关矩阵，而排序用的是绝对值，所以负边也会在“绝对值最大”的情况下被保留，这正是 4.3 要讨论的问题。

4.3 负连接、加权网络与二值化：三种常见策略

结构协变网络的偏相关矩阵几乎一定会出现负值，表示两个区域在组水平上有此消彼长的趋势（控制其他区域后）。负边的处理有三种常见策略，我分别说清楚适用场景。

策略一是最常见的“只保留正边”：把负相关视为无连接，置零后按密度阈值筛选。这个做法的优点是解释简单——“网络代表共变的区域集合”，缺点是信息有损失，且密度阈值基于绝对值选择时，强负边和强正边会互相竞争，导致正边比例被压低。

策略二是“负边单独成网”：分别用正相关和负相关构建两张网络，各自做阈值和拓扑分析。这个做法信息保全最好，但分析量和复杂性翻倍，适合探索性研究。

策略三是“保留加权符号”：不二值化，直接用带符号的加权矩阵计算全局效率等指标。缺点是有符号加权网络的图论指标解释困难，审稿时容易被追问。

我的默认做法是策略一，但在 Supplementary 里用策略二做敏感性分析。如果加权的定义是“连接强度”，那就不要在二值化时丢失权重。二值网络适合做模块检测，加权网络更适合做 hub 分析；两种网络可以同时报告，相互验证。

另一个值得注意的点：密度阈值的选择范围不是越长越好。5%–30% 是比较公认的区间，低于 5% 网络会碎成孤立节点，高于 30% 网络趋于全连接，组间差异被稀释。如果两组在任何密度下都无差异，多半是前面第 3 章的某个环节出了问题，而不是密度没调对。

5. 结构协变网络避坑指南：五个我踩过的坑

5.1 全脑几乎全连通：年龄和 ICV 混入伪协变

现象：构建出的网络矩阵相关系数普遍偏高，几乎所有 ROI 之间都有显著正相关，网络密度在 30% 阈值下依然接近满连接，且与临床分组无关。

原因：最常见的是没有回归协变量。年龄对全脑 GMV 的效应是全局性的，ICV 的效应同样如此，这两个变量会在所有 ROI 之间制造共同变异——如果两组被试的年龄或 ICV 分布不同，这些共同变异直接变成网络里的伪边。

解决：回到设计矩阵，用第 3.2 节的残差化流程，把年龄、性别、ICV 全部回归掉，再用残差构建矩阵。回归后如果网络密度显著下降，说明原来那张图里的“强连接”大部分是伪迹。

5.2 偏相关矩阵非正定：ROI 太多、被试太少

现象：调用np.linalg.inv(cov)时抛出 LinAlgError，或者偏相关矩阵出现 NaN、对角线之外出现大于 1 的数值，网络指标在相邻密度水平间剧烈跳跃。

原因：协方差矩阵不可逆，本质上是样本量不足。假设你有 148 个 ROI，却只有 50 个被试，协方差矩阵的自由度远不够，求逆结果就是数字噪音。

解决：三个选项按优先级排——减小节点数（从 148 降到 68）；改用 Ledoit-Wolf 收缩估计（见第 3.3 节）；或者直接用 Pearson 相关矩阵做网络而不是偏相关。在论文里如实报告样本量和收缩参数即可，不需要羞愧，但绝不能硬凑一个正定矩阵。

5.3 用功能连接工具直接算 SCN：输出全 1 或全 0

现象：把每个被试的 GMV 图丢进功能连接分析工具，期望它自动输出“结构网络”，得到的结果要么是单位矩阵（每个体素和自己相关），要么是两张图之间全 0。

原因：功能连接工具的设计前提是“每个被试有长时间序列”，它天然把每个体素当成一个时间点来算相关。而结构协变网络没有时间序列，只有一个静态值——工具根本不知道该算什么。

解决：结构协变网络不需要连接工具，自己构造“被试×ROI”矩阵后直接用 scipy 的pearsonr或 3.3 节的代码计算。如果实在想用现成工具，选那个工具里明确标着“group covariance”或“structural covariance”的功能，别拿普通 FC 模块硬套。

5.4 模板标签错位：左右半球和皮层下标签对不上

现象：网络矩阵中，左半球和右半球的同名区域之间的相关异常低，或者某些边和已知解剖关系完全矛盾（比如左小脑和右额叶之间出现极强相关）。

原因：分区模板的标签文件顺序与统计表的行列顺序不一致。表面分割的注释文件（annot）和皮层下分割的标签文件来自两套不同的命名体系，直接以字符串顺序合并时，很容易把第 i 行当成第 i 个 ROI，但实际对应的是另一个 ROI。

解决：合并前必须检查模板自带的 lookup table，确认每个 ROI 名称在表头中的唯一索引；并对合并后的矩阵做可视化验证——画出矩阵热图，检查对称结构的分布；再单独抽查一个解剖上就应该强相关的 ROI 对（如左右海马）验证数值合理。标签问题在神经影像里是最隐蔽的坑，因为它不报错，只给你一个错得有理有据的网络。

5.5 重测信度翻车：同一批被试两次扫描，网络对不上

现象：用同一批被试的间隔几个月两次扫描构建两套 SCN，网络拓扑指标（全局效率、模块化）差异巨大，分组差异在两次扫描间方向都反了。

原因：分割指标本身有一定噪声，头动、扫描仪状态、分割流水线版本都会引入随机变异。如果样本量不大，这些随机变异足以改变网络矩阵里中等强度的边。

解决：如果有重测数据，按第 6 章做信度筛选；没有重测数据，就把质控标准收紧——剔除分割质量差的被试，对表面指标做适度平滑（推荐 10–15mm FWHM），在报告中明确标注扫描仪型号和参数一致性。比任何统计技巧更重要的是：不要用不同扫描参数的两批数据直接拼进同一个 SCN。

6. 验证网络：信度筛选、空间置换检验与行为关联

6.1 重测信度筛选

结构协变网络构建完，第一件事不是跑模块检测，而是问一句“这条边可重复吗”。有重测数据时，我会对每条边做组内相关系数（ICC），只保留 ICC 高于 0.4 的边进入最终网络。实际操作中，ICC 低的边往往是那些平均相关强度很弱、介于阈值边缘的边，删掉它们能让结果更干净。

没有重测数据，退而求其次的做法是“分半信度”：把样本随机分成两半，分别构建网络，计算两半网络矩阵的相关，重复 100 次，取平均相关。这个相关系数如果低于 0.5，基本说明你的网络还不太稳定，优先回看分割质控和协变量回归。

6.2 空间置换检验：排除平滑和空间自相关的假象

结构协变网络的一个常见批评是：相邻脑区共享平滑效应，本来就会产生协变，网络结果可能只是空间自相关的副产品。要回应这个问题，可以用空间置换检验（spin test）：把皮层表面的形态指标图做随机空间旋转（保持空间结构、打乱位置），重复构建网络，得到“零分布”，再比较你观测到的网络指标是否显著偏离随机。

这个检验对表面数据实施比较方便，对体素数据需要用某种空间位置重排方案。不需要每个项目都做，但如果在审稿中被问到“你的网络是不是平滑导致的”，这就是标准答案。

6.3 快速验证流程：从矩阵到结论的完整链条

给新人的建议是，拿到网络矩阵后依次做三件事：第一步，检查网络密度是否在合理区间且两组差异不显著（除非研究问题本身就关心密度）；第二步，计算核心拓扑指标（全局效率、模块化、节点强度），看两组差异方向和效应量是否符合文献；第三步，把网络指标与行为/临床量表做相关，验证“网络差异有行为后果”。三步走完，结论才算闭环。

指标计算可以用 Python 的图论库，也可以用神经影像圈常用的开源图论工具包。需要注意版本兼容：有的工具包要求输入为二值网络，有的要求加权网络，传之前先看文档。我个人的习惯是保留未阈值化的加权矩阵和密度阈值化后的二值矩阵两份文件，前者用于 hub 分析，后者用于模块检测。

这个流程我跑过不止一次，最大的教训永远是：网络好不好，在分割阶段就决定了，后面的一切调参都弥补不了质控的疏忽。如果你构建的 SCN 结果让人看不懂，先别怪阈值和算法，回头去看那几个被试的分割切片图。希望这份踩坑记录能帮你少走一段弯路。

本文还有配套的精品资源，点击获取