简介这套由MATLAB编写的开源提取工具面向脑影像与神经科学研究人员用于从人类连接组项目HCP-MMP1.0分区图谱中快速提取指定脑区的体积二进制掩码。HCP-MMP1.0图谱包含精细的皮层脑区划分但直接操作图谱文件需处理复杂的体素索引工具内置脑区查找表用户只需记录目标ROI编号并调用主函数即可批量生成掩码。提取时可选择左侧、右侧或双侧区域支持将多个ROI分别保存或合并为一个大ROI并可选膨胀处理灵活适配fMRI预处理、功能连接计算、ROI统计建模等场景。资源包共7个文件核心包括主函数脚本、膨胀工具、脑区编号表、图谱压缩包和示例体积文件另附说明文档与开源许可证整包仅834KB轻量易用解压后即可使用。目前已有546人浏览学习适合需要基于HCP-MMP1.0自定义感兴趣区掩码的研究人员通过该脚本可直接得到标准MNI空间下的二值掩码省去手工绘制与坐标转换环节显著提升ROI提取效率与实验可重复性。1. 从 HCP-MMP1.0 里抠掩码fMRI ROI 分析的第一步其实最磨人每回接到 fMRI 数据处理任务第一步十有八九是从 HCP-MMP1.0 分区图谱里把目标脑区抠成二值掩码。图谱本身好拿360 个皮质区域整整齐齐可真要提取某几个区域做 ROI 分析就得写循环、对标签、处理格式最后出来的文件还不一定带正常头信息。这个开源 MATLAB 脚本把整套流程收敛成一个函数输入图谱路径和区域标签输出对应的 NIfTI 掩码文件下游不管是提取时间序列、算功能连接矩阵还是做可视化都直接复用。适合正在做 fMRI 入门研究的学生、嫌重复轮子碍事的分析人员也适合想把图谱工具移植到其他分区的开发者。2. 项目架构与提取逻辑360 个区域标签如何变成一张可用的 ROI 掩码2.1 HCP-MMP1.0 图谱的数据结构与空间约定HCP-MMP1.0 是一套基于多模态影像数据划分的人脑皮层分区图谱T1、T2 结构像和任务态、静息态功能像共同参与最终得到 360 个皮质区域左右半球各 180 个。官方发布的 volume 格式定义在标准的 MNI152 空间里通常提供 2mm 和 1mm 两个分辨率版本文件内每个体素的整数值就是区域标签背景为 0。脑影像里处理这套图谱有两种主流方式surface 版本跑在标准网格上volume 版本直接落在三维体素空间。对绝大多数做体素级 fMRI 分析的人来说volume 版是更顺手的那条路。它不需要跨网格映射直接读进来、提取、写回就能和功能像叠加。问题在于很多人拿到图谱直接load之后就开始等值比较等结果画出来才发现区域对不齐、方向翻转或者头信息丢了。要避免这个坑得先搞清楚图谱文件里到底装了什么一个三维矩阵、一套 NIfTI 头信息原点、方向、体素尺寸还有一个我们肉眼看不见的整数标签分布。理解这三件事后面的提取才有依据。2.2 文件目录结构与模块职责把这个脚本仓库解压开典型布局长这样。我拆的这份结构清爽每个文件都有明确用途fmri-extract-HCP-mask/ ├── extractHCPMask.m # 主函数单区域、多区域、组合掩码提取 ├── demo_extract_mask.m # 演示脚本自带参数可直接修改后跑通 ├── README.md # 数据下载地址和基本用法说明 ├── labels/ │ └── HCPMMP1_Label.txt # 360 个区域名与标签值的对照表 └── data/ # 图谱文件存放区体积大由用户自行下载 └── HCPMMP1_volume_2mm.nii.gz主函数只有一个文件这是刻意的设计。目标用户是脑影像研究者不是专业程序员拆成十几个.m文件反而制造认知负担。extractHCPMask.m内部自带参数校验和错误提示demo_extract_mask.m把最常用的调用方式写成注释齐全的样例照着填路径就能跑。labels目录里的对照表值得多说一句它不只是人类可读的文档脚本会读取它来做区域名到整数值的映射保证你查到的区域名和文件里的整数值对得上。data目录在仓库里是空的这很合理。HCP-MMP1.0 的体积文件在 2mm 分辨率下也有几十 MB放进仓库既臃肿又可能有版权问题。README 里通常给了官方下载指引你下载后放进data目录即可。2.3 核心处理链路读图谱、查标签、算掩码、写文件提取的逻辑链路可以稳定地分成四步每一步都有它自己的坑。第一步读图谱。niftiread返回的三维矩阵是已完成斜率和截距校正的数据也就是真实体素值。若图谱是 float 而非 int等值比较前最好先round一次避免浮点噪声把边界体素弄丢。第二步查标签。调用unique(atlas(:))把实际存在的标签值拉出来和标签对照表交叉核对。这一步看似多余实际能拦截一半以上的低级错误——我见过太多次从对照表第 n 行直接当标签值用的情况结果提取出的掩码跟预期脑区差了十万八千里。第三步算掩码。对每个目标标签做一次atlas label得到逻辑矩阵再转成 uint8内存占用比 double 小四倍写盘后文件也更小。第四步写文件。这一步最关键的是把niftiinfo读出的头信息原样回传给niftiwrite。NIfTI 头里的 qform、sform 记录着体素到现实空间的映射丢了它掩码在别的软件里打开位置就飘了。常见翻车现场就是MATLAB 里看着没问题一放进其他查看器就错位。2.4 为什么不直接用一条 shell 命令搞定有人会问很多工具里一条命令就能出掩码为什么还要用 MATLAB 脚本原因有三。第一不是所有实验室都有完整的环境而 MATLAB 在脑影像数据处理里普及度极高。第二命令行式写法没法做复杂的标签映射和校验多区域、组合掩码场景要写一串拼接可读性差。第三这个脚本的输出天然携带完整头信息并提供了前缀命名机制方便和下游 MATLAB 分析管线直接接轨。工具本身不复杂但解决的是一个高频重复需求省下的时间很可观。2.5 三种调用形态单区域、多区域、组合掩码主函数设计成向量化输入天然支持三种调用形态。单区域是最朴素的需求传一个标量标签多区域传入行向量或列向量函数内部循环生成一整套命名独立的掩码组合掩码是把多个标签做逻辑或合并成一张比如在分析左右半球对称网络时非常常用。这三种形态在 demo 脚本里都有对应样例。上手时先跑单区域确认输出文件在软件里位置正确再上多区域批量提取最后才根据自己的实验设计做组合掩码。从简到繁的节奏不容易翻车。3. 运行实操与参数拆解三步拿到第一份可校验的 NIfTI 掩码3.1 环境依赖与版本选择运行环境要求不高。MATLAB R2020a 及以上版本安装了 Image Processing Toolbox 就能直接用niftiread和niftiwrite。如果你用的 MATLAB 版本较旧、没有这两个函数也可以用脑影像社区里流行的 NIfTI 工具包替代只改读取和写回两个调用点主逻辑可以不动。操作系统方面 Windows、Linux、macOS 都跑过没有平台相关的代码。内存上不必担心。以 2mm 分辨率的 volume 图谱为例图像尺寸大约是 91×109×91体素数在百万级别uint8 读进来也就几 MB。即便用 0.7mm 的高分辨率版本也没到需要担心内存的程度——真正吃内存的是后面把掩码和四维功能像叠在一起处理时。3.2 参数对照表参数类型示例值说明atlasFile字符串data/HCPMMP1_volume_2mm.nii.gz图谱文件的绝对或相对路径labelIds数值数组[10 20] 或 10要提取的区域标签集outputDir字符串output/masks掩码输出目录不存在则自动创建prefix名值对HCP输出文件前缀默认 ROIsaveType名值对.nii输出后缀可选 .nii 或 .nii.gzprefix是个容易被忽略但很实用的参数。一个项目里如果会跑多轮提取比如健康组和患者组分别做掩码用prefix区分批次比事后改文件名靠谱得多。.nii和.nii.gz的选择看下游大多数工具都兼容老版本软件对.nii.gz支持一般通常选.nii更稳。3.3 三步跑通改路径、执行 demo、核对输出第一步准备好图谱文件放到data目录里。第二步打开demo_extract_mask.m把atlasFile、labelIds、outputDir改成自己的实际值。下面是 demo 里最核心的三行几乎不需要动其他东西% demo_extract_mask.m 中的关键调用 atlasFile fullfile(data, HCPMMP1_volume_2mm.nii.gz); % 换成你的图谱路径 labelIds [4 10 20]; % 示例提取三个区域 outputDir fullfile(output, masks); % 输出目录自动创建 extractHCPMask(atlasFile, labelIds, outputDir, prefix, demo);参数含义atlasFile指向图谱文件labelIds里每个数字都是图谱中实际存在的标签值建议先跑一次unique确认outputDir不存在时函数内部会mkdir不会因为目录缺失中断prefix, demo让生成的三个文件分别叫demo_label_4.nii、demo_label_10.nii、demo_label_20.nii。第三步在 MATLAB 里切到项目根目录运行demo_extract_mask。终端会逐标签打印处理进度像这样Processing label 4 ... done (12345 voxels) Processing label 10 ... done (23456 voxels) Processing label 20 ... done (34567 voxels) 3 mask(s) written to output/masks.看到每行voxels数是合理范围几千到几万不等取决于区域大小基本可以认为提取正确。3.4 输出校验从文件大小到头信息文件生成后别急着进下游。我用三条命令做快速体检每次都执行几乎不花时间% 体检 1文件是否存在且非空 d dir(output/masks/demo_label_10.nii); fprintf(File size: %.2f KB\n, d.bytes / 1024); % 体检 2掩码非零体素数是否和终端输出一致 data niftiread(output/masks/demo_label_10.nii); fprintf(Nonzero voxels: %d\n, sum(data(:) 0)); % 体检 3头信息里的体素尺寸是否保留 info niftiinfo(output/masks/demo_label_10.nii); fprintf(Voxel size: %.2f x %.2f x %.2f mm\n, ... info.PixelDimensions(1), info.PixelDimensions(2), info.PixelDimensions(3));第一项检查落盘成功第二项验证标签值换算正确第三项确认空间元数据没丢。脑影像数据最怕静默错误这三行体检能拦住绝大多数问题。如果PixelDimensions输出和原始图谱一致就可以放心把掩码喂给下游流程了。3.5 命令行批量运行不打开图形界面也能跑如果要把提取流程放进批处理脚本里MATLAB 的-batch模式是好选择。一条命令就能跑无 GUI适合在服务器环境或者多被试批量分析时使用matlab -batch extractHCPMask(data/HCPMMP1_volume_2mm.nii.gz, [4 10 20], output/masks, prefix, demo)注意-batch模式下 MATLAB 不会进入交互环境脚本路径相关的问题更容易暴露所以命令行调用时建议传入绝对路径避免工作目录不清楚导致文件找不到。这个模式对跑多被试批量数据很有用每个被试一个输出目录循环调用即可。4. 核心代码拆解从图谱读取到掩码写回的完整链路4.1 主函数入口与参数解析% extractHCPMask.m % 从 HCP-MMP1.0 分区图谱中提取任意区域的三维二值掩码 % 参数说明 % atlasFile : 图谱 NIfTI 文件支持 .nii 和 .nii.gz % labelIds : 单个或多个区域标签例如 10 或 [10 20] % outputDir : 输出目录不存在时自动创建 % 可选参数prefix 前缀saveType 输出格式 function extractHCPMask(atlasFile, labelIds, outputDir, opts) arguments atlasFile string {mustBeFile} labelIds {mustBeNumeric} outputDir string opts.prefix string ROI opts.saveType string .nii end % 读取图谱及完整头信息 vol niftiread(atlasFile); info niftiinfo(atlasFile); % 校验图谱必须是三维或四维且存在有效标签 assert(ndims(vol) 3 ndims(vol) 4, ... 图谱维度异常请检查输入文件);arguments块是 MATLAB 较新版本里的正式参数用法好处是类型错误在进函数体之前就会被拦下来。mustBeFile校验路径存在mustBeNumeric防止有人把字符串当标签传入。第三件重要的事info在读完文件后立刻接上取出来为后面写掩码保留原始头信息。很多人在这一步偷懒直接用niftiwrite(mask, outFile)不传头信息结果掩码空间参考全部丢失。4.2 标签校验与掩码生成% 列出图谱中实际存在的标签值避免标签文件与图谱不一致 actualLabels unique(vol(:)); % 检查请求的标签是否都真实存在 missing setdiff(labelIds, actualLabels); if ~isempty(missing) warning(以下标签在图谱中不存在: %s, mat2str(missing(:))); end % 逐标签生成二值掩码并写出 for k 1:numel(labelIds) label labelIds(k); fprintf(Processing label %d ... , label); % 等值比较生成逻辑掩码再转为 uint8 节省内存 mask (vol label); mask uint8(mask); % 组合掩码文件路径 outFile fullfile(outputDir, ... sprintf(%s_label_%d%s, opts.prefix, label, opts.saveType)); % 写回前显式声明 uint8 类型保持头部与数据一致 info.Datatype uint8; niftiwrite(mask, outFile, info); fprintf(done\n); end end逻辑说明vol label生成了和vol尺寸相同的逻辑矩阵目标体素为true其余false转成uint8后数值变成 1 和 0这正是“二值掩码”的定义。这里特别注意setdiff的处理标签不存在时我只给 warning 不中断因为多区域批量提取场景下漏掉一两个区域不该让整个任务失败。那行assert(ndims(vol) 4)的校验针对的是四维图。若有人不小心传入四维功能像等值比较仍能运行但会按时间帧全部匹配掩码语义就不对了所以提前拦截。写回前info.Datatype uint8是容易被忽略的细节默认继承原图谱的 datatype 可能会把掩码写成 float下游查看器兼容性就差。4.3 多标签组合一张掩码覆盖多个脑区% combineHCPMasks.m % 把一组标签合并成一张组合掩码适用于左右半球对称 ROI 分析 function combineHCPMasks(atlasFile, labelIds, outFile) vol niftiread(atlasFile); info niftiinfo(atlasFile); % 初始化全零掩码 mask zeros(size(vol), uint8); % 对每个标签做逻辑或写操作保留一次图谱读取 for k 1:numel(labelIds) mask(vol labelIds(k)) 1; end % 体积可能含第四维写文件前压缩到三维 mask squeeze(mask(:, :, :, 1)); info.Datatype uint8; niftiwrite(mask, outFile, info); fprintf(Combined mask written: %s\n, outFile); end这里的关键不是等值比较本身而是“一次读取多次扫图”。vol只被加载一次内部循环在内存中反复跑到同一份数据时间开销远小于每轮重新读取磁盘。mask(vol labelIds(k)) 1用逻辑索引直接落点不需要临时变量三个标签合并成一张掩码时这段代码的可读性就体现出来了。4.4 等值比较与 ismember 的选择对单个标签vol label是最快的写法。如果标签是一组ismember(vol, labelIds)能一次性生成多标签掩码但代价是内存占用瞬时升高——MATLAB 会先构造匹配结果矩阵高分辨率图谱下多出几倍内存。脚本里默认单标签循环是因为 2mm 图谱的循环开销极低而我见过不少人在 1mm 图谱上硬用ismember一次提取五十个区域内存直接打到十几 GB。折中方案就是 4.3 的组合函数外层循环、内层逻辑索引兼顾速度和内存。5. 避坑与常见问题标签偏移、空间重采样和文件格式的实战记录脚本能跑通和跑出正确结果之间隔着不少边界情况。下面几条是我拆这个项目时遇到过、或者身边同事反复踩过的坑每条都按现象、原因、解决的顺序写希望对你有实际帮助。5.1 标签偏移从对照表查到的是名字不是标签值现象按照HCPMMP1_Label.txt第 n 行的区域名提取生成的掩码覆盖区域和预期明显不符有的甚至是另一个半球的位置。原因对照表行号不一定等于 volume 文件里的标签值有些重打包版本还做过整体重编号行号和标签值的偏差可能不是简单的 1而是乱序。解决提取前先在 MATLAB 跑unique(vol(:))把实际标签打印出来和对照表核实再设定labelIds。我用这套流程处理过多个来源的图谱每次都能发现差异已成固定动作。如果你的区域名是字符串建议先写一个小脚本把对照表里的名字映射到标签值再传给主函数双保险。5.2 空间不匹配图谱空间和功能像空间不是一回事现象生成的掩码在功能像上叠加时错位或者提取的时间序列噪声巨大相邻 voxel 的值忽高忽低。原因HCP-MMP1.0 volume 图谱定义在 MNI152 标准空间与未经配准的 native 空间功能像不同坐标系形状差异大。解决先做空间标准化把功能像配准到图谱空间后再用掩码或者反过来把图谱重采样到功能像分辨率并改成最近邻插值重采样后对数据先 round 再取整避免中间灰值干扰。我一般偏向先把功能像统一到图谱空间因为这样图谱始终保持整数标签不动提取逻辑最干净。5.3.nii.gz读取失败压缩格式不是每个 MATLAB 版本都受支持现象niftiread报错说文件扩展名不支持或者读出来全是零数组。原因较老的 MATLAB 版本对.nii.gz支持不完整某些第三方下载的文件虽然后缀是.nii.gz内部 header 字段却不标准。解决先gunzip转成.nii再喂给脚本如果仍报错用 hex 编辑器看文件开头字节确认是 NIfTI-1 还是 NIfTI-2 格式对应换工具包。注意转换后记得把原始.nii.gz保留一份别一手滑删掉原始数据。5.4 输出文件在外部软件里错位或全黑现象掩码在查看器里位置漂移或者在软件里整体黑屏看不见任何体素。原因写文件时没复用info头信息或数据类型被默认写成了 double部分查看器对 double 掩码支持得并不好。解决niftiwrite之前先uint8(mask)写回时传info并显式设置info.Datatype uint8打开文件后检查信息面板里的 voxel size 是否正常。若手头有原始图谱和掩码做一次体素间距离对齐测试几秒钟就能判断头信息对不对。5.5 批量提取 360 个区域时越跑越慢现象循环 360 个标签前面飞快后段越跑越慢甚至假死。原因没有把vol保持在内存每次迭代重新niftiread或输出目录写在共享盘上大量小文件的 IO 成为瓶颈。解决把读取提到循环外4.3 的代码就是这么组织的输出先写本地临时目录全部成功后movefile到共享盘再不行就tic/toc定位是哪一步慢。脑影像数据的小文件操作比想象中伤性能不要忽视 IO。6. 二次开发从掩码提取到 ROI 时间序列与多图谱扩展6.1 组合掩码与左右半球对称提取利用 HCP-MMP1.0 左右半球标签的偏移规律可以一行代码生成左右对称的组合 ROI% 提取左右半球对称组合掩码 leftLabels [10 20]; rightLabels leftLabels 180; % 左侧 1-180右侧 181-360 allLabels [leftLabels, rightLabels]; combineHCPMasks(data/HCPMMP1_volume_2mm.nii.gz, allLabels, ... output/combined_symmetric.nii);参数说明leftLabels是左半球两个区域的标签180得到右侧对应区域。这个技巧在构建双侧网络 ROI 时很常用五秒钟就能完成过去要手写两轮提取的活。6.2 从掩码到时间序列下游联动掩码最终要服务于分析。下面是提取 ROI 平均时间序列的常用写法逐帧计算内存友好% 读取四维功能像并提取 ROI 平均时间序列 roiMask niftiread(output/demo_label_10.nii) 0; funcData niftiread(func_preprocessed.nii); % 四维 [x y z t] nFrames size(funcData, 4); ts zeros(nFrames, 1); for t 1:nFrames frame funcData(:, :, :, t); ts(t) mean(frame(roiMask)); end逐帧循环对时间序列长度几百帧的情况开销可接受但避免了把整个四维数组按掩码展开成二维矩阵后者在高分辨率功能像上非常吃内存。拿到ts之后就可以做功能连接、组间比较、或者作为特征输入后续模型了。6.3 适配其他图谱这个提取函数本身对图谱内容无感知核心就是“读图谱、查标签、算掩码、写文件”四步。换成其他分区图谱时只需要换一套标签对照表、确认目标图谱的空间定义主函数一行不用改% 适配其他分区图谱标签定义在 1-170 区间 extractHCPMask(data/other_atlas_1mm.nii.gz, [31 32], output/other, ... prefix, OTHER_ROI);从那以后我每次拿到一套图谱文件都会先做一次标签直方图和头信息体检确认空间参考、体素尺寸、datatype 三项无误再写任何提取逻辑。这套习惯帮我少踩了不少坑也让掩码在下游流程里几乎一次通过。希望帮到你。本文还有配套的精品资源点击获取
