简介：针对小鼠宽场光学成像（WOI）数据处理与统计分析的开源MATLAB工具箱，由神经影像团队开发并公开验证。资源面向神经科学、功能磁共振成像及光学成像领域的研究人员、博士生和技术开发者，旨在弥补WOI数据缺乏统一分析框架的空白。文章详细介绍了工具箱的设计原理、具体函数模块与统计方法，整合了多组WOI和fMRI技术，支持细胞特异性钙指示剂信号处理及类血氧动力学数据分析，并能在光栓塞中风模型中可靠检测感觉缺陷、绘制电爪刺激激活区。全文还阐述了去噪、平滑等预处理流程以及多重比较校正策略，有助于提高信噪比和检测能力。包体为单个PDF文件，大小1.43MB，内容结构完整，配有公开数据集说明，便于读者按模块复现与应用。已有127人学习，适合从事小鼠脑功能网络、疾病模型及神经影像分析的读者参考。
1. 宽场光学成像缺的不是采集设备，而是一套标准化数据处理工具箱
宽场光学成像（Wide-field Optical Imaging, WOI）能在同一套系统上同时记录小鼠皮层血流动力学信号和细胞特异性钙信号，空间分辨率达到约 100 μm/像素、时间分辨率 16.8 Hz，这正好填补了 fMRI 与双光子显微镜之间的尺度空白。但真正让实验室卡住的不是采集，而是数据处理：数据拉出来是 pixels×pixels×frames 的三维堆栈，从去噪、掩模、滤波到统计检验，每一步都缺少像 fMRI 社区那样成熟的工具链。这篇发表在 Neurophotonics 的论文，把 fMRI 里经过十几年验证的聚类分析加随机场理论翻译到 WOI 像素空间，做成了一个开源 MATLAB 工具箱，并用光栓中风模型和电爪刺激实验验证了整套流程。对正在跑 WOI 数据、或者打算从 fMRI 转做光学成像的研究组来说，这套工具箱可以直接作为数据处理框架的起点，省掉从零造轮子的过程。
2. 从原始数据到大脑掩模：load_data 与解剖标志点的正确打开方式
2.1 数据格式与载入：先想清楚帧结构再动手
工具箱的第一个入口是 START 文件夹下的 load_data.m。它做的事情很朴素——把成像系统导出的数据读入 MATLAB，并抽出一帧作为后续掩模绘制的参考图。但这一步里藏着一个容易忽略的前提：数据必须是 pixels×pixels×frames 的三维堆栈，而不是把多波长通道混进去的四维数组。
% load_data.m 核心逻辑：读入三维数据堆栈并抽取参考帧 [filename, filepath] = uigetfile('*.mat', '选择WOI数据文件'); data = importdata(fullfile(filepath, filename)); % 数据格式要求：像素 x 像素 x 帧数，多波长需提前拆分 assert(ndims(data) == 3, '数据必须是三维堆栈，请检查是否混入波长通道'); frame5 = data(:, :, 5); % 论文取第5帧作为掩模绘制的参考帧第 5 帧这个参数不是随手填的。成像开始的前几帧往往伴随 LED 切换或荧光漂白的不稳定，取中段帧能获得相对平稳的皮层强度分布。如果你自己的系统稳定帧偏后，就把 frame5 换成晚几帧的索引，这个变量就是给你留的调整口。另一个容易踩的点是多波长数据：论文系统有 470/530/590/625 nm 四个 LED 依次激发，如果直接把原始导出文件丢给 load_data.m，三维断言会直接报错，必须先按波长拆开、每个波长单独存一个三维堆栈。
2.2 掩模怎么画才不返工：roipoly 与大脑边界的判断
拿到参考帧后，make_mask_landmarks.m 会调用 MATLAB 内置的 roipoly.m，让你沿着皮层轮廓点击一圈，生成一个二值掩模——大脑区域为 1，非大脑区域为 0。这个掩模后面会乘到每一帧上，把颅骨、缝合线、窗片边缘全部滤掉，只保留有效像素参与统计。
% make_mask_landmarks.m 核心逻辑：交互式绘制大脑掩模 figure; imagesc(frame5); axis image; colormap gray; % roipoly 是MATLAB内置函数，点击闭合区域后返回二值掩模 isbrain = roipoly(); % 掩模会保存为变量 isbrain，后续处理乘到整段数据上绘制掩模看起来简单，但实际操作有讲究。参考帧的亮度在整个 FOV 内并不均匀——光线从液导进入皮层时中心亮、边缘暗，olfactory bulb 和 superior colliculus 区域强度低，容易让人误以为那不是大脑。我的做法是先把显示窗口的对比度拉满，确认脑膜血管走向，再顺着皮层凸面的亮带勾边。如果你画的掩模明显偏离皮层凸面曲率、或者把两侧颞肌的反射区框进去了，说明参考帧显示窗口没调好。返回值 isbrain 是决定后续所有统计量的门卫，画歪一点，功能连接图就多一片鬼影。
2.3 两个解剖标志点决定整个分析的几何基准
掩模画完之后，make_mask_landmarks.m 里还有一步关键交互：程序会弹框让你点击前囟缝（anterior suture）和 lambda 两个解剖标志点。这两个点不是摆设——它们一方面用于生成种子点原型，另一方面为后续仿射变换提供坐标基准。
从这两个坐标出发，工具箱会自动在皮层表面铺一套种子区域，对应运动皮层、体感皮层、视觉皮层等经典分区，并弹出预览图。如果种子原型左右不对称、或者没有覆盖整个背侧皮层表面，选择“否”可以重新点标志点，直到对称性满意为止。
% MakeSeedsMouseSpace.m 核心逻辑：根据标志点生成种子原型 [x_suture, y_suture] = ginput(1); % 点击前囟缝 [x_lambda, y_lambda] = ginput(1); % 点击lambda % 根据两点连线与距离，外推标准种子模板到当前图像空间 % 然后显示预览，由用户判断是否对称、是否覆盖全皮层这一步是整个管线里最依赖操作者判断的环节。两个点的连线方向决定种子模板的旋转角度，两点距离决定尺度缩放。同样是 C57BL/6 小鼠，不同个体的颅骨标志点间距会有零点几毫米的差异，如果你只凭感觉点、不确认种子预览就点“是”，后面把所有小鼠对齐到公共空间时会发现形变量特别大。论文的流程里允许点“否”重新来过，这正是为了强制你检查这一步的几何质量。建议每只小鼠、每次手术后在同样放大率下采集参考帧，这样标志点的相对位置差异会被控制到最小。
3. 系统相关与系统无关处理：光谱校正、信号提取和滤波的先后顺序
3.1 proc1_sys_dep.m 到底在干什么：暗场、光谱响应与 LED 波长
处理管线的第一步是 proc1_sys_dep.m，系统相关步骤。做宽场成像的人对这一步又爱又恨——它处理的是硬件给你留下的坑，换了成像系统必须重新校准。
第一步是暗场扣除。采集一套 Dark.tif（关闭激发光时传感器记录的底噪），从数据中减去，避免 CMOS 暗电流被当成真实的荧光信号。第二步是光谱校正：论文系统用了 Mightex LED 和 Semrock 二色镜，每个波长的光在传播路径上都有不同的衰减，需要用 prahl.txt 这类光谱文件来归一化不同波长的响应差异。这一步的意义在于，不让光学系统本身的透过率曲线扭曲你真正关心的生理信号。
% proc1_sys_dep.m 核心逻辑：暗场扣除与光谱归一化 dark = imread('Dark.tif'); % 每套系统需单独采集暗场 data_dark = data - dark; % 扣除底噪 % 光谱归一化：根据prahl.txt中的透过率曲线逐波长校正 % 470nm激发峰处的响应权重最高，其余波长按比例衰减 data_corrected = data_dark ./ spectral_response;这段话值得刻在屏幕上：Dark.tif 是每套系统专属的。论文用的是 Zyla 5.5 sCMOS，你如果换了相机、换了镜头、甚至换了液导光纤，暗场和光谱响应曲线都要重新采集。我见过有人直接拿论文配套的 Dark.tif 用在自己的系统上，结果暗电流不匹配，后期功能连接的噪声基底整体抬高，这是典型的“参数跟着论文走”翻车案例。
3.2 Proc2.m 的信号提取逻辑：为什么系统无关
处理完系统相关步骤，接下来是 Proc2.m，系统无关处理。这个模块处理的是数据本身的结构性问题，包括数据裁剪、坏帧剔除和信号平滑。
信号平滑这一步需要特别说明。fMRI 里常做空间平滑，因为体素之间需要借邻域信息提升信噪比；WOI 数据因为组织散射，本身相邻像素之间就存在有效 FWHM 达数个像素宽的光学模糊，再做空间平滑的尺度就要克制。论文的 Proc2.m 更强调时间方向的预处理：检查数据是否有异常跳变帧，必要时去掉头尾不稳的帧段，保证进入统计的是一段平稳信号。此外，这一步还会把掩模乘上去，让后续的计算只在 isbrain 为 1 的像素里进行。
% Proc2.m 核心逻辑：系统无关的信号整理 data_masked = data_corrected .* isbrain; % 掩模应用到每一帧 % 检查坏帧：全局信号突变超过3个标准差则标记 bad_idx = find(abs(zscore(mean(data_masked, [1 2]))) > 3); data_clean = data_masked; data_clean(:, :, bad_idx) = []; % 剔除坏帧 whos data_clean; % 确认清理后的维度坏帧阈值 3 个标准差是通用做法，但要注意一个坑：电爪刺激实验里，刺激期间的全局信号变化可能本身就超过 3 个标准差。如果你把刺激期当成坏帧剔除，激活响应就没了。我一般先看全局均值曲线，确认哪些跳变是刺激锁时的、哪些是突发的伪迹，再决定阈值。
3.3 滤波频段不是随手填的：钙信号与血流信号的频率分工
滤波是整个管线里物理意义最明确的一步，也是最容易用错的一步。论文里给出了两组滤波参数——钙信号做 0.4 到 4.0 Hz 的 Butterworth 带通，对应 delta 频段；血红蛋白信号做 0.009 到 0.08 Hz 的带通，对应 infraslow 频段。
这两组参数不是随便定的。GCaMP 钙信号的生物学节律集中在低频但高于血流信号，0.4 Hz 以下主要是血流动力学串扰，4.0 Hz 以上是生理噪声（呼吸、心跳在麻醉小鼠身上常常落在 5 Hz 以上）；而血红蛋白信号走的是和 BOLD 一致的 infraslow 通路，0.009 Hz 以下接近直流漂移，0.08 Hz 以上则开始混入钙信号的泄漏成分。
% OIS_GCaMP_Filter.m 核心逻辑：频段分离滤波 fs = 16.8; % 采样率，论文系统每波长16.8 Hz % 钙信号：delta频段 0.4-4.0 Hz calcium_filt = bandpass(data_clean, [0.4 4.0], fs); % 血红蛋白信号：infraslow频段 0.009-0.08 Hz hemo_filt = bandpass(data_clean, [0.009 0.08], fs);麻醉小鼠心率约 300-600 次/分，对应 5-10 Hz；呼吸频率约 90 次/分，对应 1.5 Hz。1.5 Hz 其实落在钙信号频段内，所以如果你的数据呼吸伪迹明显，就需要额外加一个陷波或更高的滤波阶数。Butterworth 滤波器的阶数默认设置如果不够，通带边缘会有明显抖动，看频谱图能快速判断：如果通带内有孤立尖峰，那就是生理噪声没滤干净，别急着跑统计，先回来处理这个尖峰。
这里必须强调滤波方向。数据是时间序列，MATLAB 的 bandpass 默认零相位滤波，双向通过，这没问题。但如果你用 filtfilt 之前的旧习惯 split 处理，相位偏移会直接扭曲 FC 计算里像素间相关的时间对齐关系——相关性算出来的不是神经同步而是滤波假象。工具箱里用的是带通滤波器函数，默认零相位设计，不需要额外担心这一点。
4. 功能连接怎么算才有意义：种子点、双侧连接和运动质控
4.1 种子点 FC 的像素级相关计算
预处理完的数据，下一步就是算功能连接（FC）。工具箱里的功能连接分析分成种子点 FC 和双侧 FC 两套，各有各的应用场景。
种子点 FC 的做法是：以第 2 章生成的种子区域为参考，提取每个种子区域的平均时间序列，然后和整个 FOV 内每一个像素的时间序列做皮尔逊相关。得到的结果是一张相关图，每个像素的值代表该像素与种子区域的同步程度。落回到小鼠模型上：体感皮层种子应该和同侧的感知运动网络同步，而不是和视觉皮层强相关。
% seed_FC.m 核心逻辑：种子点与全像素的相关分析 seed_ts = squeeze(mean(data_filt(:, :, seed_mask > 0), 1)); % 提取种子区平均时序 for i = 1:numel(isbrain) if ~isbrain(i), continue; end [r, p] = corrcoef(squeeze(data_filt(i)), seed_ts); fc_map(i) = r(1, 2); % 存储像素与种子的相关系数 end这里有个统计口径要注意。相邻像素因为组织散射的存在，有效空间分辨率大约覆盖多个像素，像素间天然有空间相关性。这意味着 FC 图里很大一片显著区域可能只是光学模糊的功劳，而不是生物学上的功能连接。所以单独看一张 FC 图意义不大，关键是通过后续的聚类统计来校正——这在第 5 章展开。
4.2 双侧 FC 图：左右对应是天然的参照
双侧 FC 是另一种分析方式，计算左右半球的对称区域之间的连接强度。这招在小鼠脑成像里特别实用：小鼠皮层两半球在解剖上高度对称，正常情况下左右体感皮层之间应该有稳定的同步。论文的光栓模型做的是左侧体感前爪区，3 天后左侧病灶区域与右侧的同源区域连接强度显著下降——这就是双侧 FC 的直接读出。
工具包的实现方式是：以中线为镜面，把一侧的种子模板翻转映射到另一侧，逐区域计算两侧时间序列的相关。这个分析对仿射变换的依赖更高——如果两只小鼠之间没有对齐到公共空间，镜像映射就完全错位，算出来的双侧 FC 会低得离谱。
4.3 check_mvmt.m 这条质控线不能跳过
运动伪迹是活体光学成像的宿敌。小鼠虽然做了头固定，但呼吸、咀嚼、甚至血管搏动都会让皮层图像出现亚毫米级位移。工具箱里的 check_mvmt.m 用时间导数的全局方差来评估帧间运动量：计算每一帧与下一帧的差，然后看整个 FOV 内的方差。
% check_mvmt.m 核心逻辑：估计帧到帧的全局运动响应 temporal_diff = diff(data_filt, 1, 3); % 时间导数 global_var = squeeze(mean(var(temporal_diff, [], 1:2))); plot(global_var); % 应该是一条平稳的基线，有明显尖峰说明有运动理想情况下 global_var 是一条平稳线，如果某个时间点出现明显尖峰，那一段数据就要被标记甚至剔除。关键经验是：越到数据后期，运动伪迹越容易悄悄变大——小鼠在固定装置里待久了，肌肉紧张度变化导致缓慢漂移，这种慢漂移在全局方差曲线上不一定表现为尖锐峰，而是一段缓慢抬升的平台。做 FC 分析前必须看这条曲线，而不是只看有没有大尖峰。
5. 避坑：多重比较校正的“翻译”过程常见翻车点
5.1 为什么 Bonferroni 在 WOI 里寸步难行
做像素级统计的人都知道多重比较校正的痛。一张 156×156 的 FC 图，除掉掩模外区域，FOV 内剩下约一万多个像素。如果对每个像素做 t 检验，就是一万多次假设检验。Bonferroni 校正把所有检验当成相互独立的，把显著性阈值除以检验次数——在 WOI 数据上这条路走不通：相邻像素在物理上被组织散射耦合，天然不独立，Bonferroni 会把一堆真实激活区全部扼杀。
论文翻译了 fMRI 社区的标准做法：先设定一个像素级阈值（如 p < 0.001），找出超过阈值的像素，再计算这些像素组成的聚类大小，通过与随机场理论（RFT）导出的空分布对比，判断这个聚类是否显著。本质上是把“单个像素显著与否”的问题升级为“邻接像素组成的团块是否足够大”——大团比孤立像素更有统计说服力。
% cluster_threshold.m 核心逻辑：根据RFT估计聚类大小阈值 % 关键输入：FWHM（图像有效分辨率）、FOV内像素数、像素级阈值 fwhm_pixels = 2.5; % 组织散射造成的有效FWHM，单位：像素 p_thresh = 0.001; % 像素级阈值 cluster_thr = compute_cluster_threshold(fwhm_pixels, p_thresh, num_pixels); % 输出：超过 cluster_thr 像素的数目的聚类才被判定为显著FWHM 的估计值直接决定了阈值大小。论文分析中使用了从光传输仿真得到的有效点扩散函数宽度，如果你自己系统的 FWHM 和文献差异很大，比如换成高倍物镜或更薄的颅窗，就必须重新估算 FWHM——这个参数错 10%，聚类阈值可能错 30%。
5.2 五条血泪经验
现象一：同样的数据用 Bonferroni 校正，FC 激活区几乎全部消失；换成聚类阈值校正，体感激活区域干干净净显示出来。
原因：相邻像素高度相关，Bonferroni 的独立性假设完全不成立，它把效应量分散到一万多个相关检验上，真实信号被淹没了。
解决：采用论文翻译的聚类数阈值。用 FWHM 捕获像素间残余相关性，再用团块大小而非单像素 p 值做决策。fMRI 社区用这套办法积累了近二十年，WOI 数据的空间平滑特性其实比 fMRI 更适配它。
现象二：每只小鼠的种子模板位置看起来都对，但群体统计时功能连接图一片糊。
原因：不同小鼠解剖标志点的相对位置有差异，直接在原始空间做平均，相当于把不同坐标系的功能图硬套在一起。
解决：跑 Affine.m 做仿射变换，把所有小鼠映射到一个公共模板空间。这个步骤必须放在 FC 计算之前，因为种子模板的位置随小鼠变化，会影响所有下游统计。
现象三：滤波后的钙信号里出现周期性噪声条纹，抑制不住。
原因：带通滤波后的信号混入了呼吸或心跳的谐波成分，尤其是 1.5 Hz 的呼吸基频和它 3 Hz 的二次谐波，后者的耦合峰正好穿过钙信号频带。
解决：先看功率谱再定滤波器参数。如果谱上有孤立窄峰，重建一个级联的陷波滤除，或者调整带通的频段边缘——从 0.4-4.0 Hz 改成 0.4-3.5 Hz 可以把二次谐波甩出通带，代价是损失一点高频波动信息，但在实际荧光测量里不相干噪声的抑制更值得。
现象四：电爪刺激实验里，对照组一半的“激活区”出现在对侧镜像位置。
原因：没有扣除刺激带来的系统性伪迹，例如 LED 激发强度的波动与刺激器同步发生，被相关分析解读成“功能连接”。
解决：检查 check_mvmt.m 的全局方差曲线，看刺激期间是否有全局跳变；同时看刺激对照条件是否有类似的“激活区”，如果有，说明是系统性的伪迹，得回到 proc1_sys_dep.m 查暗场和光谱校正是否到位。
现象五：同一个数据集跑两遍，聚类显著的激活区大小差很多。
原因：掩模绘制和标志点选择两次不一致，导致参与统计的像素数变了，FWHM 估计也漂了。
解决：把掩模保存下来，整个研究用同一份掩模和同一组标志点坐标；如果需要对不同批次实验，建议只取一只代表动物的标志点，其余全部通过仿射变换对齐到这个基准。
6. 进阶：仿射对齐、批处理与结果验证的实操技巧
6.1 Affine.m 把不同小鼠对齐到公共坐标空间
跨个体统计分析的前提是空间对齐。论文的管线里 Affine.m 负责把每只小鼠的图像仿射变换到共同的 atlas 空间，算法以第 2 章选定的前囟缝和 lambda 两点作为锚点，计算刚性（旋转、平移、缩放）变换矩阵，把数据重切片到公共网格。
% Affine.m 核心逻辑：基于标志点计算仿射变换矩阵 A = estimate_affine(src_pts, ref_pts); % src_pts为默认参数 data_aligned = imwarp(data_filt, A, 'linear'); % 线性插值重切片做这一步时有两个细节值得注意。第一，插值方法选择线性而不是最近邻——最近邻插值在重切片时会引入锯齿状像素级伪影，直接污染后续的聚类大小估计；第二，变换前后 FOV 的大小会变化，需要裁剪到公共区域，不然不同小鼠的有效像素范围不一致，统计效能在边缘区域会明显下降。
6.2 复现论文实验：光栓中风模型的缺陷检测
论文用两个实验验证了整套框架。第一个是光栓中风模型——在左侧体感前爪皮层用 532 nm 激光配合 Rose Bengal 诱导光栓，3 天后用工具箱检测到该区域的功能连接缺陷；第二个是电爪刺激激活图，用同样的流程清晰分离出刺激对应的皮层激活区。
如果你想复现这两个实验，最值得借鉴的是论文 Day 0 和 Day 3 的自身对照设计。基线（Day 0）数据是同一只小鼠的“健康参照”，Day 3 的病变数据直接对比，这样个体差异这个最大的噪声源天然被抵消。群体的 FC 变化分析时，每个像素位置算 Day 3 减去 Day 0 的差值，再做一次 t 检验和聚类校正，得到的就是病变直接导致的功能变化图。这个差值策略比直接组间比较省力得多，也能显著减少需要的动物数量。
6.3 批量处理的习惯
工具箱里单只小鼠的处理流程是清晰隔离的，但实验常常要跑十几二十只鼠。我会把整批处理的组织方式提到最高优先级——为每只小鼠建一个独立文件夹，记录掩模文件、标志点坐标和滤波参数快照。每只动物跑完基础处理后，强制生成一张质控图（参考帧、掩模轮廓、种子模板、全局方差曲线四联图），归档进同一个文件夹。
从那以后我每次批量处理都强制走一遍这个流程：先做一只小鼠全流程确认参数，再顺手归档质控图，然后才把同一套流程原封不动地应用到剩下的小鼠。这种死板的做法看起来慢，但避免了“全部跑完后发现其中几天采集的 Dark.tif 换了”这样的灾难性返工。希望帮到你——至少这套工具箱让 WOI 的分析不至于卡在你这边。
本文还有配套的精品资源，点击获取
