1. 为什么还写这套SPM12的批处理脚本1.1 它解决什么问题MATLAB配合SPM12做fMRI预处理算是神经影像老牌组合了。这两年虽然fMRIPrep、Nipype这类工具越来越流行但很多课题组的老数据、旧脚本、以及正在跑的纵向研究仍然跑在SPM12这套流程上。我的情况跟很多同行类似手头有一大批.nii.gz格式的fMRI数据需要跑完从时间层校正到空间平滑这一整套标准预处理而实验室服务器上又只有MATLAB和SPM12没有配Python环境更没有fMRIPrep容器。这时候一套自己写的批量预处理脚本就成了刚需。这套脚本解决的核心问题有两个。第一把SPM12图形界面里点点点的操作变成命令行批量执行几十个被试一次性跑完不用守在电脑前手动点上千次按钮。第二处理最常见的压缩格式.nii.gzSPM12自身识别这种压缩格式偶尔会出问题所以脚本里要先做解压和整理把这些“脏活”一并自动化。它不是什么高深的东西但在只有MATLAB的老环境下这套办法是最省事、最可靠的一条路。1.2 适合谁来看这篇东西如果你正准备用SPM12做一批fMRI数据的预处理尤其是手头有旧版本SPM12、数据是压缩格式、或者需要批量处理多个被试这篇内容应该能帮你省下不少折腾的时间。我下面写的脚本思路和代码不是拿过来就能复制的成品而是一个可以照着改、照着拼的框架。你需要根据自己的数据实际情况调整TR、扫描层数、参考图像、分割模板路径等参数但整体流程和坑点我都替你踩过一遍了。2. 数据准备与预处理流程设计2.1 数据目录先整理对后面才省心拿到一批.nii.gz格式的数据第一件事不是写预处理脚本而是先把数据目录整理成统一的结构。我最开始图省事数据直接从服务器上拷下来一股脑丢进一个文件夹结果后面的批量循环脚本写得一脸懵因为不同被试的文件命名规则不统一。后来老老实实按下面的结构重排了一遍一切顺畅多了。project_dir/ raw/ sub-01/ func/ sub-01_task-rest_run-1_bold.nii.gz sub-01_task-rest_run-2_bold.nii.gz anat/ sub-01_T1w.nii.gz sub-02/ func/ ... anat/ ... preprocess/ sub-01/ sub-02/这一步花不了多少时间但能帮你少掉至少三成调试脚本的精力。凡是涉及批处理的任务目录结构就是地基地基不稳后面全白搭。实际整理的时候我会先写一个极简的MATLAB脚本读入原始文件夹下的子目录列表然后用movefile把每个被试的functional和anatomical文件分别挪到对应位置。这一步不需要什么高级技巧核心逻辑就是两步读目录、挪文件。实际上我更推荐在处理之前用Python的os.rename来做但如果你只有MATLAB用movefile也一样只是记得在循环里用fullfile拼路径。2.2 预处理流程的六个核心环节SPM12的fMRI预处理一般围绕下面六个环节展开。我用的这套旧版SPM12流程顺序是固定的每个环节之间的依赖关系也很明确。时间层校正Slice TimingfMRI扫描是分层采集的不同层面的获取时间点不一样所以需要校正到同一个参考时间点。这一步需要明确的参数有两个一个是TR重复时间即扫描一个完整脑体积所需的时间另一个是扫描层数以及各层之间的采集顺序。如果你的数据是隔层采集这里填错了后面的时间序列分析基本就废了。头动校正Realign被试在扫描仪里不可能完全不动头动会造成不同时间点上同一体素采集到的组织不一致。Realign的核心操作是把所有时间点上的图像配准到第一个时间点或者平均图像上输出六个头动参数三个平移、三个旋转。这六个参数后面做回归分析时也是重要协变量所以一定要保存。配准Coregister这一步是把功能像和结构像对齐让fMRI数据能叠加到高分辨率的T1像空间里。核心是选择参考图像通常是T1像和源图像通常是平均功能像其他文件跟着源图像一起变换。这里最常见的坑是忘记把source文件选对导致功能像和结构像没对齐。分割Segment把T1结构像分割成灰质、白质、脑脊液等组织类别。SPM12的分割过程还会同时计算一个从个体空间到标准空间的变形场deformation field这个文件是下一步归一化必需的。归一化Normalize利用上一步生成的变形场把每个被试的功能像和结构像都配准到标准空间MNI空间。这一步的关键是选择“写出的文件”write而不是重新估计参数否则会重复计算变形场浪费时间。空间平滑Smooth对归一化后的功能像做高斯核平滑提高信噪比。最常见的参数是FWHM半高全宽一般设置为[6 6 6]但具体要看你的体素大小和分析目的后面我会详细说怎么选。一段流程走完每个人脑都落在同一个标准空间里不同被试、不同组之间的数据就有了可比性。研究套路深不深另说这一步的标准程度直接影响后面统计结果的可靠性。3. 批量预处理脚本的核心实现3.1 脚本整体骨架先看一下我用的脚本整体结构。这个脚本不是一次写完的中间改了很多版目前这套是在MATLAB R2016b和SPM12旧版r7219环境下跑的。如果你用的是更新版本个别函数接口可能有差异但整体逻辑完全可以移植。% 主脚本: spm12_batch_preprocess.m % 功能: 批量fMRI预处理 (.nii.gz) for SPM12旧版 % 环境: MATLAB R2016b, SPM12 r7219 clear; clc; % 设置SPM路径 spm_path /path/to/spm12; addpath(spm_path); spm(defaults, FMRI); % 设置fMRI默认参数 spm_jobman(initcfg); % 初始化job管理器 % -------------------------- % 1. 配置全局参数 % -------------------------- TR 2.0; % 重复时间单位秒 num_slices 32; % 扫描层数 slice_order interleaved; % 隔层采集 ref_slice 16; % 参考层通常取中间层 smooth_kernel [6 6 6]; % 平滑核大小 data_root /data/project; output_root /data/project/preprocess; % 获取所有被试列表 subject_list dir(fullfile(data_root, raw, sub-*)); n_subjects length(subject_list); % -------------------------- % 2. 逐被试处理 % -------------------------- for subj 1:n_subjects fprintf(Processing %s (%d/%d)\n, subject_list(subj).name, subj, n_subjects); subj_dir fullfile(data_root, raw, subject_list(subj).name); out_dir fullfile(output_root, subject_list(subj).name); if ~exist(out_dir, dir); mkdir(out_dir); end % 2.1 解压 .nii.gz gunzip_files(subj_dir); % 2.2 生成各步骤的matlabbatch run_preprocessing_pipeline(subj_dir, out_dir, TR, num_slices, slice_order, ref_slice, smooth_kernel); end fprintf(All done.\n);这段主脚本的核心逻辑就是一个循环加上两个自定义函数gunzip_files负责解压run_preprocessing_pipeline负责生成并运行各步骤的SPM批处理任务。后面我会把这两个函数的关键代码展开。3.2 解压.nii.gz这一步为什么不能省SPM12在读取.nii.gz压缩格式时不同版本的支持情况不一样。我用的这套旧版经常出现“cant open file”或者读到一半报错的情况所以稳妥起见第一步先把所有.nii.gz解压成.nii再让SPM去读。解压函数没什么技术含量但有几个细节值得注意。第一解压前先检查目标文件是否已经存在避免重复解压浪费磁盘IO。第二解压后要删除原来的.nii.gz文件吗我的建议是别删留着原始备份万一处理过程出问题还能重新来过。第三解压后的文件命名要和原文件严格一致只去掉.gz后缀不额外改动文件名。function gunzip_files(subj_dir) % 递归找出所有.nii.gz文件 file_list dir(fullfile(subj_dir, **, *.nii.gz)); for i 1:length(file_list) src_file fullfile(file_list(i).folder, file_list(i).name); [path, name, ext] fileparts(src_file); % ext .gz [~, base_name] fileparts(name); % 去掉.nii后缀 target_file fullfile(path, [base_name, .nii]); % 如果已经解压过跳过 if exist(target_file, file) fprintf( [skip] %s already exists\n, target_file); else fprintf( [gunzip] %s\n, src_file); gunzip(src_file); end end end用fileparts处理文件名时要注意name变量拿到的其实是“xxx.nii”还需要再套一次fileparts才能拿到“xxx”这个纯文件名。这个细节我一开始没注意解压出来的文件全是乱名后来用[~, base_name] fileparts(name)处理了一下才解决问题。3.3 每个预处理步骤的matlabbatch怎么写SPM的批处理本质上是构造一个叫matlabbatch的MATLAB结构体然后调用spm_jobman(run, matlabbatch)执行。这个结构体长得跟SPM图形界面上点出来的配置是一一对应的所以最保险的写法是先在图形界面上手动配置一次然后保存batch文件直接看SPM生成了什么结构照着改就行。下面是run_preprocessing_pipeline函数的核心代码我按步骤拆开讲。步骤一Slice Timingmatlabbatch{1}.spm.temporal.st.scans { {func_files} }; % 所有功能像文件每一卷作为一个cell matlabbatch{1}.spm.temporal.st.nslices num_slices; matlabbatch{1}.spm.temporal.st.tr TR; matlabbatch{1}.spm.temporal.st.ta TR - TR/num_slices; matlabbatch{1}.spm.temporal.st.so slice_order_vector; % 实际采集顺序例如1:2:32, 2:2:32 matlabbatch{1}.spm.temporal.st.refslice ref_slice; matlabbatch{1}.spm.temporal.st.prefix a;这里有个参数值得一提tatemporal acquisition表示一次完整体积采集时间内实际用于采集的时间长度。公式是ta TR - TR/num_slices这个值SPM不会自动帮你算必须手动填。忘了填或者填错的后果就是时间层校正的结果会整体偏移看起来处理完了实际上校正了个寂寞。slice_order_vector需要根据实际的采集方式生成。对于隔层采集interleaved常见序列是奇数层先采、偶数层后采对应向量应该是[1:2:num_slices, 2:2:num_slices]。如果是连续采集sequential就是1:num_slices。这个信息可以从扫描仪的参数里查到或者用dcm2niix转换数据的时候输出信息里一般会带着。步骤二Realignmatlabbatch{2}.spm.spatial.realign.estwrite.data { {func_files} }; matlabbatch{2}.spm.spatial.realign.estwrite.eoptions.rtm 1; % 重对齐到平均图像 matlabbatch{2}.spm.spatial.realign.estwrite.eoptions.interp 4; % 4阶B样条插值 matlabbatch{2}.spm.spatial.realign.estwrite.eoptions.wrap [0 0 0]; matlabbatch{2}.spm.spatial.realign.estwrite.eoptions.quality 0.9; matlabbatch{2}.spm.spatial.realign.estwrite.roptions.which [2 1]; % 输出所有图像平均图像 matlabbatch{2}.spm.spatial.realign.estwrite.roptions.interp 4; matlabbatch{2}.spm.spatial.realign.estwrite.roptions.wrap [0 0 0];这里的rtmregister to mean建议设成1也就是把每一时间点的图像对齐到所有时间点的平均图像上比对齐到第一幅图像更稳定。roptions.which我写的是[2 1]意思是写出的文件包括重新采样后的所有图像2和平均图像1。平均图像后面做Coregister的时候要用到。步骤三Coregistermatlabbatch{3}.spm.spatial.coreg.estwrite.ref { {t1_file} }; matlabbatch{3}.spm.spatial.coreg.estwrite.source { {mean_func_file} }; matlabbatch{3}.spm.spatial.coreg.estwrite.other {func_files_resliced}; % 所有功能像 matlabbatch{3}.spm.spatial.coreg.estwrite.eoptions.cost_fun nmi; matlabbatch{3}.spm.spatial.coreg.estwrite.eoptions.interp 4; matlabbatch{3}.spm.spatial.coreg.estwrite.roptions.interp 4; matlabbatch{3}.spm.spatial.coreg.estwrite.roptions.which [0 1];Cost function我习惯用nmi互信息它对多模态图像T1和fMRI的配准效果比ncc归一化互相关稳定。other字段里不仅要填所有功能像还要包含前面Slice Timing和Realign输出的文件。Coregister计算的时候会以source为基准估计变换参数然后把source和other里的所有文件一起按这个变换写入新的空间。步骤四Segmentmatlabbatch{4}.spm.spatial.preproc.channel.vols {t1_file}; matlabbatch{4}.spm.spatial.preproc.channel.biasreg 0.001; matlabbatch{4}.spm.spatial.preproc.channel.biasfwhm 60; matlabbatch{4}.spm.spatial.preproc.tissue(1).tpm {fullfile(spm_path, tpm, TPM.nii,1)}; matlabbatch{4}.spm.spatial.preproc.tissue(1).ngaus 2; % ... 组织类别2到6依此类推 matlabbatch{4}.spm.spatial.preproc.warp.mrf 1; matlabbatch{4}.spm.spatial.preproc.warp.cleanup 1; matlabbatch{4}.spm.spatial.preproc.warp.bb [-78 -112 -70; 78 76 85]; matlabbatch{4}.spm.spatial.preproc.warp.vox [1.5 1.5 1.5]; matlabbatch{4}.spm.spatial.preproc.warp.write [1 1]; % 写变形场和反变形场Segment步骤里最需要注意的是TPM组织概率图的路径。旧版SPM12的TPM在spm12/tpm/TPM.nii每个组织类别是一卷所以引用时要带上,1、,2这样的卷序号。很多报错都出在这个地方——路径不对或者卷标不对整条处理链直接中断。warp.write [1 1]表示同时写出正向和反向变形场。正向是把个体空间映射到MNI空间的变形场后面Normalize用反向是把MNI空间映射回个体空间的变形场通常留作备用。步骤五Normalizematlabbatch{5}.spm.spatial.normalise.write.subj.def {deformation_field_file}; % 上一步生成的y_*.nii matlabbatch{5}.spm.spatial.normalise.write.subj.resample {func_files_to_normalize}; matlabbatch{5}.spm.spatial.normalise.write.woptions.bb [-78 -112 -70; 78 76 85]; matlabbatch{5}.spm.spatial.normalise.write.woptions.vox [3 3 3]; % 功能像归一化到3mm体素 matlabbatch{5}.spm.spatial.normalise.write.woptions.interp 4; matlabbatch{5}.spm.spatial.normalise.write.woptions.prefix w;Normalize步骤里同样要注意的是bbbounding box和vox体素大小这两个参数。功能像归一化通常用[3 3 3]或[2 2 2]结构像用[1 1 1]。它们只需要在MATLAB的job配置里设定一次SPM会自动应用到你指定的所有文件上。步骤六Smoothmatlabbatch{6}.spm.spatial.smooth.data {normalized_func_files}; matlabbatch{6}.spm.spatial.smooth.fwhm smooth_kernel; matlabbatch{6}.spm.spatial.smooth.dtype 0; matlabbatch{6}.spm.spatial.smooth.im 0; matlabbatch{6}.spm.spatial.smooth.prefix s;平滑的fwhm是核心参数。它设置的是高斯核的半高全宽单位是毫米。常见做法是设置成跟体素大小差不多的值比如体素3mm时用[6 6 6]体素2mm时用[4 4 4]。这个参数不宜过大过大会把一些本来就有意义的精细结构也抹掉也不宜过小过小起不到平滑的作用。smooth.dtype 0表示输出数据类型跟输入一致一般就是float32够用。im 0表示不生成隐含掩膜这个默认参数不用改。所有步骤的matlabbatch组装完成后用多行spm_jobman(run, matlabbatch)依次执行每一步。我建议分步执行而不是把六个job一股脑塞给spm_jobman(run, matlabbatch)一次跑完原因有两点一是可以在每两步之间检查输出文件二是如果某一步参数写错你马上就能定位到具体是哪一步出的问题不至于整个batch全部失败后从头排查。3.4 参数计算和选型实战以TR、FWHM为例上面提到的好几个参数看起来就是填数字但实际设置时是需要根据数据算的。我以TR和FWHM为例具体说说怎么算、怎么选。TR的计算很简单一般来说就是从扫描仪或者dicom头文件里读出来的值。假设你的扫描参数是TR2000ms那就是TR2.0。填入SPM时单位是秒不是毫秒。如果TR信息丢失了可以用一个小技巧估算看单个.nii文件的时间维长度和总扫描时长。比如一个volume有300个时间点总扫描时长是10分钟那TR≈600/3002秒。这只是估算如果有原DICOM文件还是以DICOM头里的Repetition Time字段为准。FWHM的选择就比较讲究了。它跟你的体素大小直接相关同时跟你的分析目标也有关系。做组水平的统计分析时一般建议FWHM取两倍体素大小左右。例如体素3mm取6mm合适体素2mm取4mm。但如果你做的是高空间分辨率的单个被试分析或者用MVPA这类保持空间细节的方法平滑核就应该小一些甚至不做平滑。一个实用的判断标准是看你的fMRI数据本身的空间分辨率平滑核大小不应超过感兴趣区域尺寸的一半。从参数到代码这些数值都是直接填进上面的matlabbatch结构里的。我习惯把参数统一定义在主脚本开头用一个变量名指代这样后面要改只需要改一处不用翻遍脚本找数字。3.5 只处理部分被试怎么加断点续跑跑批量预处理经常遇到的情况是跑到第10个被试时服务器崩了或者某个被试数据质量差导致spm_jobman报错中断。如果脚本没有断点机制前面的9个被试就白跑了浪费时间。我在主循环里加了一个简单的断点判断为每个被试建立一个“完成标记文件”。for subj 1:n_subjects subj_name subject_list(subj).name; done_flag fullfile(output_root, subj_name, finished.txt); if exist(done_flag, file) fprintf([skip] %s already processed.\n, subj_name); continue; end % ... 处理该被试 fid fopen(done_flag, w); fprintf(fid, completed %s\n, datestr(now)); fclose(fid); end这样做的逻辑很简单每次跑完一个被试就写一个标记文件下次再跑脚本时自动跳过已经完成的被试。同样的思路可以精细化到每个步骤比如给Slice Timing、Realign都单独写标记这样中断后可以从具体步骤重新续跑。我在实际项目中通常会用“每跑完一个步骤手动把对应job从循环里注释掉”这种最土的办法但效果也很实在。4. 常见问题与排查实录4.1 SPM报错类型和解决方法我把这几年来用这套脚本时遇到的高频报错整理成了一个表方便你对照排查。现象可能原因解决方法Cant open file xxxx.nii.nii.gz未解压或路径中有中文/空格先跑解压步骤路径全英文无空格Error using spm_run_normalise变形场文件没生成或者y_*.nii文件名匹配失败检查Segment步骤是否真的完成确认deformation_field_file变量路径Undefined function spmSPM12路径没加进MATLAB搜索路径检查addpath(spm_path)是否执行成功确认SPM是完整解压Out of memory功能像文件太大或同时加载了太多文件减少一次加载的文件数分批处理run增加虚拟内存Job execution failedSPM版本不一致导致matlabbatch字段不兼容大部分旧版脚本的字段多一个write子结构对照图形界面重新保存新版本batchImage dimensions do not match不同run或不同被试图像的维数不一致检查原始数据确认时间层数和矩阵大小一致后才能一起跑Realign第一个问题最常见。很多初学者从网盘或课题组服务器拷到数据压缩包后直接塞给SPM处理结果SPM读不了.nii.gz。用我这个脚本里的gunzip_files函数先解压问题就解决了一大半。第二个问题也很有代表性。SPM12的Segment步骤会生成一个y_*.nii变形场文件但旧版SPM12生成的变形场文件名前缀不同如果你用的是更新版本这个文件名前缀可能变成y_或者别的导致下一步Normalize找不到文件。我的建议是在脚本里加一行检查def_file spm_select(FPList, out_dir, ^y_.*\.nii$); if isempty(def_file) error(Deformation field not found. Check Segment step.); end用spm_select的FPList模式可以按正则表达式从指定文件夹里筛选文件比手动拼接路径靠谱得多。类似的在后面步骤里筛选w*、sw*文件也都可以用这个方式。4.2 头动参数过大的数据要不要剔除Realign跑完之后每个被试都会生成一个rp_*.txt文件里面是六列头动参数。我通常会在预处理完成后写一个小段代码批量读取所有被试的头动参数计算最大平移和最大旋转然后标记那些超过阈值的被试。% 读取头动参数 rp_file spm_select(FPList, subj_dir, ^rp_.*\.txt$); motion load(rp_file); max_trans max(abs(motion(:, 1:3))); % 单位mm max_rot max(abs(motion(:, 4:6))); % 单位rad % 旋转角度转换为度 max_rot_deg max_rot * 180 / pi; fprintf(%s: max_trans%.2fmm, max_rot%.2f度\n, subj_name, max(max_trans), max(max_rot_deg));我的经验阈值是最大平移超过3mm或者最大旋转超过3度的被试就需要谨慎处理。这不意味着一定要剔除但至少要在后续统计分析时作为协变量纳入或者检查一下是不是存在明显的尖峰。如果数据本身质量不好头动超过阈值可以试着先用ART工具或者tsdiffana之类的质量检查工具跑一遍看有没有异常时间点再决定是剔除整个被试还是去掉某些时间点。这已经属于预处理后面的质量检查环节了但建议在跑完整套流程后立即做不要等统计分析时才发现数据没法用。4.3 旧版SPM12的几个特殊问题“旧版”这两个字在SPM12的语境下意味着你的代码和某些新版本的行为可能存在差异。我用的是r7219这一版遇到的主要差异有三个。第一matlabbatch结构在旧版接口更繁琐。比如新版SPM12的Realign步骤有时简化了字段旧版还保留着estwrite.eoptions和estwrite.roptions这些层级你只要照着图形界面保存的batch结构去写就行不要凭新版教程照搬。第二TPM模板文件的卷标路径格式。旧版SPM12里TPM模板需要写成fullfile(spm_path, tpm, TPM.nii,1)这样的格式有些版本能够自动识别不加卷标的写法但旧版不行必须显式带上,1。这个细节导致过很多次“Index exceeds matrix dimensions”之类的报错排查起来很隐蔽。第三旧版对“脏数据”更敏感。比如说.nii文件头缺少某些必要字段新版会补默认值旧版直接报错。遇到这种情况可以用dcm2niix重新转换原始DICOM文件确保NIfTI头信息完整再跑预处理。4.4 并行加速matlabpool和parfor能省多少时间批量处理的痛点在时间。一个被试跑完六步机械硬盘环境下大概需要15分钟到半小时。30个被试就是8到15小时基本过夜。如果服务器是多核CPU可以尝试用MATLAB并行计算来加速。最直接的办法是给主循环换一个parfor。但这里有个前提每个被试的SPM任务之间不能共享任何工作区变量matlabbatch结构体必须完全独立构建。我的做法是把run_preprocessing_pipeline这个函数写成无副作用的纯函数只输入路径和参数输出状态信息然后就可以放心用parfor。parfor subj_idx 1:n_subjects % 注意不要用subject_list(subj_idx).name这种动态变量名最好先把属性提取出来 subj_name subject_names{subj_idx}; fprintf(Processing %s (%d/%d)\n, subj_name, subj_idx, n_subjects); subj_dir fullfile(data_root, raw, subj_name); out_dir fullfile(output_root, subj_name); if ~exist(out_dir, dir); mkdir(out_dir); end gunzip_files(subj_dir); run_preprocessing_pipeline(subj_dir, out_dir, TR, num_slices, slice_order, ref_slice, smooth_kernel); endparfor的坑也不少。最典型的是每个worker都要初始化SPM——SPM本身不是并行安全的如果两个worker同时初始化SPM经常会互相冲突。解决办法是在run_preprocessing_pipeline函数内部调用一次spm(defaults, FMRI)和spm_jobman(initcfg)让每个worker独立初始化。我实测过的场景是服务器上8个物理核心并行跑30个被试从过夜跑缩短到大约3小时。提速很明显但代价是内存占用大尤其Smooth和Segment步骤。如果服务器内存不够大建议并行数控制在4核左右2个被试一起跑比强行开8个核然后把内存撑爆划算。5. 脚本的后续扩展方向如果你只是要跑完标准六步预处理上面的内容完全够用了。但实际项目中经常会碰到需要额外处理的情况。我把自己实验后觉得值得扩展的几个方向列在下面。加一个ART或DVARS检查环节预处理完成后用Artifact Detection Tools或者自己写一个小脚本计算DVARSframe-wise displacement的导数识别异常时间点。这可以作为后续回归分析的头动协变量来源也可以用来筛选被试。加入去噪步骤SPM12旧版自带的DARTEL工具可以改进分割和归一化精度。如果你的数据样本量不大、时间预算充足建议在Segment和Normalize之间插入DARTEL创建模板的步骤。这样做的主要收益是组水平配准精度更高但流程会拉长很多。串联conn工具箱做功能连接预处理如果你的课题是做功能连接那么预处理之后通常还有一步去噪和带通滤波。CONN工具箱可以复用SPM12的处理结果直接输入sw*.nii文件即可。我自己做静息态功能连接项目时就是把这套预处理脚本的输出作为后续CONN分析的输入衔接得非常顺利。生成分析报告给每个被试生成一张预处理质量报告内容包含头动曲线、配准前后覆盖叠加图、归一化前后对比图。我后期整理文章和对外汇报时这些图省了我不大事。用SPM自带的spm_check_registration配合MATLAB的print函数就能批量导出。这些扩展方向不需要重新设计主流程都是在现有脚本框架上打补丁加模块。预处理管线这种东西只要骨架搭得稳往上面挂新组件很容易。我个人在跑完这么多批数据后最大的体会是批量预处理的核心不在于某一步参数调得多完美而在于每一步的输出都能被严格检查、每一次异常都能被快速定位。这套脚本虽然“旧”但因为每一步之间的文件传递都是显式的、可控的所以它反而比很多黑盒工具更容易排查问题。如果你正准备上手自己的fMRI数据建议先把这套流程完整跑通一批再根据结果去调整参数细节比一开始就追求“完美预处理”要靠谱得多。