动态被试间功能连接揭示由连续任务或交流范式驱动的脑网络瞬时配置

摘要

该研究方法的目的是确定在呈现时间锁定刺激的功能磁共振成像记录过程中，功能连接在何时（时间维度）以及在哪些脑区之间（空间维度）发生显著的重组变化。

摘要

基于任务的功能性磁共振成像在理解大脑对各种刺激的反应方面具有巨大潜力；然而，这种研究通常忽略了功能处理的动态特性，其分析结果往往混杂了任务驱动效应与大脑活动固有的自发波动。本文介绍了一种可突破上述局限性的新方法学流程：采用滑动窗分析方案，能够追踪功能连接随时间的变化，结合被试间的相关性测量，该方法可分离出纯粹与刺激相关的效应。得益于严格的阈值化处理过程，可以提取并分析被试间功能相关性的显著变化。

在一组接受自然主义视听刺激的健康受试者上，我们通过将揭示出的功能重组与电影中的特定线索相关联，展示了该方法的实用性。我们展示了如何通过本方法捕捉脑活动的时间动态特征（即特定连接的演变过程），或聚焦于关键时间点上的空间活动模式。我们提供了一个完全公开可用的完整分析流程，并逐步说明了其使用方法及关键参数的影响。

引言

功能磁共振成像（fMRI）已成为一种首选工具，用于无创监测由外部刺激引起的脑活动变化。更具体而言，人们日益关注对不同脑区激活时间序列之间统计依赖性的理解，这种依赖性被称为功能连接性（FC）¹，通常以皮尔逊相关系数进行计算。大量研究已表明，全脑范围内的功能相互作用会随着所执行任务的不同而发生动态重组^2,3,4。

为了超越这一初步表征，已有两个分析方向被分别推进：一方面，研究发现，在特定脑区中由时间锁定的刺激所引发的反应，在不同受试者之间表现出强烈的关联性⁵。对这种被试间相关性（ISC）进行量化，显示出在深化认知^6,7,8,9和脑部疾病^10,11理解方面的潜力。此外，这种跨被试的相关性分析方法还被进一步扩展至对跨脑区同步性的评估¹²，即后来被称为被试间功能相关性（ISFC）的方法¹³。

另一方面，功能连接（FC）重构的动态特性开始受到越来越多的关注（参见Hutchison 等¹⁴、Preti、Bolton和Van De Ville¹⁵、Gonzales-Castillo和Bandettini¹⁶关于静息态和任务态下此问题的最新综述）。特别是，通过在逐渐移动的时间子窗口上进行连续的相关性测量，可以追踪全脑功能连接随时间的变化^17,18，从而在行为任务背景下提供更深入的见解^19,20。

本文提出了一种将上述两种方法相结合的方法学框架。具体而言，我们采用滑动窗口的方式计算ISFC，以追踪受试者在经历时间锁定的自然范式时，跨脑区同步性的动态演变过程。由于该方法具有跨被试的特性，分析主要聚焦于刺激驱动的效应，而被试间不相关的自发性fMRI变化则被显著抑制。这一点至关重要，因为静息态与任务诱发的脑活动模式正逐渐被认识到具有截然不同的特性^21,22。

至于该方法的动态组分，它能够对任务刺激进行更全面且准确的表征，尤其是在探究一种自然范式时，该范式中多种线索（听觉、视觉、社会性等）随时间组合呈现。此外，由于对显著动态波动的合理统计评估一直存在激烈争论^23,24，我们的方法通过将显著的ISFC变化与适当的零假设数据进行比较，特别关注并妥善处理了分析中的这一方面。

我们以一组健康受试者为研究对象，这些受试者接受了视听电影刺激，我们展示了可以从局部电影子区间中准确提取出时间与空间上的ISFC变化特征。在此过程中，我们还描述了用户需选择的主要分析参数的影响。本研究所呈现的结果基于先前已发表数据的一部分^25,26。

方案

本方案已获得当地伦理委员会批准（Biomedical Inserm 365 方案 C08-39）。

1. 成像前准备

1. 招募研究对象群体，并获得所有受试者的书面知情同意。需获得当地伦理委员会的批准。
2. 选择一种可对所有受试者以时间锁定方式实施的研究范式。
   注意： 本研究中采用了一部面向青少年的视听科学纪录片（https://miplab.epfl.ch/index.php/miplife/research/supplement-asd-study）。

2. 成像

1. 针对分析中考虑的每位受试者，至少进行一次功能成像扫描，使被扫描的志愿者接受与时间锁定的兴趣范式。
   1. 使用3特斯拉磁共振成像（MRI）扫描仪，通过平面回波成像序列获取横断面图像。
   2. 采用以下成像参数：体素大小 = 3 mm × 3 mm × 3 mm，重复时间（TR）= 2 s，回波时间 = 50 ms，视野 = 192，40层切片。
      注意：在可行范围内建议采用更短的TR值。该方案也可使用更小的视野（例如，当分析仅限于特定脑亚结构时），这将有助于实现更高的时间分辨率（更低的TR）或空间上更精确的分析。
   3. 在刺激呈现前后保留数秒的记录时间（≥ 2个TR）。
2. 至少进行一次独立的功能成像扫描，使被扫描的志愿者在扫描仪中处于静息状态，闭眼，并被指示不要入睡。
   注意：将刺激相关和静息态数据采集分开，可避免不同条件之间可能存在的相互影响（例如，事先观看电影可能会对后续采集的静息态数据留下持续痕迹）²⁷。如果不需要进行上述额外的静息态数据采集，流程中可采用另一种计算方法（尽管该方法更易检测出假阳性；见讨论部分），即用从范式相关信号计算出的替代时间序列来代替该数据（见步骤5.1.2）。
3. 进行结构成像。
   1. 使用3特斯拉磁共振成像（MRI）扫描仪，以及T1加权磁化准备快速梯度回波序列。
   2. 采用以下成像参数：体素大小 = 1 mm × 1 mm × 1 mm，视野 = 256，176层切片。

3. 数据与软件准备

1. 针对每个待分析的会话，请确保以下数据文件的存在：
   1. 一组功能磁共振成像（fMRI）数据，以独立的三维 NIFTI 或 HDR/IMG 文件形式存在，且具有统一的编号方案（例如，“fMRI_0001”、“fMRI_0002”等）。
   2. 一幅T1加权结构磁共振成像（MRI）图像，格式为NIFTI或HDR/IMG。
   3. 一个位于蒙特利尔神经学研究所（Montreal Neurological Institute, MNI）空间的兴趣图谱，格式为NIFTI。
      注：代表性受试者（“S17”）所需的输入文件示例及完整分析流程代码可在 https://c4science.ch/source/Intersubj_pipeline.git  获取
2. 下载最新版本的公开可用的 Freesurfer 软件²⁸（https://surfer.nmr.mgh.harvard.edu/fswiki/DownloadAndInstall）。
3. 从 https://www.fil.ion.ucl.ac.uk/spm/software/spm12/ 下载最新版本的公开可用的统计参数映射（Statistical Parametric Mapping, SPM）MATLAB 工具箱。
4. 打开 MATLAB（版本 2017a 或更新版本），定位到新下载的“freesurfer”和“spm12”文件夹。对每个文件夹，右键单击并选择 添加到路径 > 所选文件夹和子文件夹 选项。

4. 数据预处理

1. 在 MATLAB 终端中输入 spm 以启动 SPM12 主菜单，然后点击 fMRI 按钮，进入专用于 fMRI 数据的预处理选项。对每个 fMRI 会话分别执行以下步骤以进行预处理。
   1. 点击 Realign (Est & Res)，在新打开的批处理编辑器（Batch Editor）窗口中，双击 Data > Session。在新弹出的 Session 窗口中，选择所有需要处理的功能像图像。然后点击 Done 按钮，接着点击批处理编辑器窗口中的 Run Batch 图标（绿色三角形）。等待重对齐步骤完成，MATLAB 终端窗口中将显示完成提示。
   2. 点击 Coregister (Est & Res)，在新打开的批处理编辑器窗口中，双击 Reference image。在新弹出的 Reference Image 窗口中，选择下一步生成的平均功能像体积文件（以 "mean" 开头），然后点击 Done 按钮。接着双击 Source image，在新弹出的 Source Image 窗口中选择 T1 图像。点击 Done 按钮，然后点击批处理编辑器窗口中的 Run Batch 图标（绿色三角形）。等待配准步骤完成，MATLAB 终端窗口中将显示完成提示。
      注意：此步骤将覆盖原始 T1 图像，使更新后的图像与功能像数据处于相同的空间中。
   3. 点击 Segment，在新打开的批处理编辑器窗口中，双击 Volumes。在新弹出的 Volumes 窗口中选择 T1 图像，然后点击 Done 按钮。接着在批处理编辑器窗口中，双击 Deformation Fields 并选择 Inverse 选项。点击 Run Batch 图标（绿色三角形），等待分割步骤完成，MATLAB 终端窗口中将显示完成提示。
2. 在 MATLAB 终端中输入 JOVE_GUI1，以打开第一个预处理图形用户界面窗口。对每个 fMRI 会话执行以下步骤以进行分析。
   1. 点击 Enter fMRI data，并选择在步骤 4.1.1 中生成的所有重对齐后的功能像文件（以 "r" 开头）。对于 IMG/HDR 文件格式，需同时选择 IMG 和 HDR 文件。
   2. 在专用的可编辑文本框中输入数据的 TR 值（单位为秒）。
   3. 点击 Enter T1 data，并选择在步骤 4.1.3 中生成的三个概率性组织类型图像（分别以 "c1"、"c2" 和 "c3" 开头）。
   4. 点击 Enter motion file，并选择在步骤 4.1.1 中生成的、对应会话的运动参数文本文件（以 "rp" 开头）。
   5. 选择所需的预处理类型，即是否对数据进行去趋势处理（分别通过设置专用单选按钮为开启或关闭状态），以及应去除哪些协变量（从专用列表中选择适当选项）。
      注意：该回归步骤借鉴自 DPARSF 工具箱中的一个函数²⁹。对每位受试者的白质和脑脊液信号，取其在 DPARSF 概率性组织模板中信号值大于 0.99 的体素进行平均。在我们的分析中，对数据进行了去趋势处理，并去除了白质/脑脊液的时间序列信号，以及常数项、线性项和二次项趋势。
   6. 点击 Preprocess 以开始预处理，并等待结果在窗口中显示。用户可通过修改选项后再次点击 Preprocess 按钮，以不同方式重新进行预处理。
      注意：灰质图像的绘制方式借鉴了 Power 等人提出的表示方法³⁰。
   7. 如需保存输出结果以供后续步骤使用，请点击 Save 按钮。如需清空窗口内容，请点击 Clear 按钮。

补充图1：首次预处理图形用户界面窗口的示例截图。 根据所选预处理选项处理后的灰质体素的逐体素时间进程（右上图），以及可用于预处理的协变量（从上至下：脑脊液/白质平均时间进程、平移运动参数和旋转运动参数）。请点击此处查看该图的放大版本。

3. 在 MATLAB 终端中输入 JOVE_GUI2 以打开第二个预处理图形用户界面窗口。对每个待分析的 fMRI 会话执行以下步骤。
   1. 点击 选择数据，选择在步骤 4.2.7 中保存的数据文件（名为“ISFC_VX.mat”）。
   2. 点击 选择运动，选择步骤 4.1.1 中生成的、以“rp”为前缀的包含运动参数的文本文件。
   3. 点击 选择图谱，选择用于分割的图谱所对应的NIFTI文件。
   4. 点击 选择逆向变换，选择第 4.1.3 步中生成的、代表从 MNI 空间到原生空间变形场的 NIFTI 文件，该文件前缀为 "iy"。
   5. 点击 选择 fMRI 体数据，并选择任意一个fMRI数据卷。
      注意： 此步骤用于获取功能数据的头部信息，因此所选的具体体积并不重要。
   6. 在专用的可编辑文本窗口中输入数据的重复时间（TR，单位为秒）。
   7. 输入与清洗相关的信息：要执行的清洗类型（即, 在“Scrubbing type”列表中选择在标记帧前后需剔除的帧数，以及帧间位移阈值（Power标准）³¹在“剔除阈值”（Scrubbing threshold）可编辑文本窗口中输入应剔除fMRI数据体的阈值（单位：mm）。
      注意： 对经过清理的数据点进行三次样条插值，以邻近样本的估计值替代。在我们的分析中，我们在标记的体积数据之后清理了一个时间帧，并采用0.5 mm的阈值进行清理。
   8. 输入用于ISFC计算的滑动窗口大小W（以TR为单位）（见步骤5）。
      注意： 该信息将允许通过源自DPARSF工具箱的函数对时间过程进行滤波处理²⁹，当 f = 1/W Hz 时³²在我们的分析中，我们采用 W = 10 TR 作为折中值，以捕捉动态波动，同时保留足够的样本用于稳健估计。
   9. 点击 绘图 单击按钮以显示去伪影和滤波步骤前（上方图）及步骤后（下方图）的指示性图谱化时间进程。通过视觉检查确认，在选定的预处理步骤之后，这些输出信号未包含明显的伪影成分。
   10. 为保存后续步骤的输出结果，请在专用的可编辑文本窗口中输入保存名称，然后点击 保存 按钮。要清除窗口内容，请单击 清晰 按钮

补充图2：第二次预处理图形用户界面窗口的示例截图。 完成图谱映射后、在（上图）和经过（下图）根据所选参数进行数据擦除（scrubbing）与滤波处理后的区域时间进程。每条曲线代表从所有可用区域中随机选取的一个区域的时间进程。请点击此处查看该图的放大版本。

5. 滑动窗口ISFC计算

1. 在 MATLAB 终端中输入 JOVE_GUI3 以打开第一个与 ISFC 相关的图形用户界面窗口。对每种类型的已获取 fMRI 会话片段（刺激相关片段、刺激相关会话的静息态片段以及纯静息态片段）分别执行以下步骤。
   1. 点击 加载数据，并选择通过步骤 4.3 创建的所有适当数据文件。
   2. 选择所选会话片段是否应进行相位随机化。
      注意：若无可用的静息态记录，相位随机化可作为从刺激相关信号生成零假设数据的替代方法。
   3. 在专用的可编辑文本窗口中输入数据的 TR（以秒为单位）。
   4. 在专用的可编辑文本窗口中输入用于分析的滑动窗参数：用于计算连接性测量的窗口大小（以 TR 为单位），以及连续窗口应移动的步长大小（以 TR 为单位）。
      注意：在我们的分析中，使用了 10 TR 的窗口大小和 1 TR 的步长大小。
   5. 修改“会话类型”表格，以指明所加载的会话片段中哪些是在相同实验条件下采集的。使用从 1 开始递增的整数标记不同类型的片段（例如，在某次记录中刺激是第一次还是第二次呈现）。如果仅采集了一种类型的会话片段，则保持该表格不变。
      注意：在本研究中，“会话”可指结合电影/静息态的记录（如 图 1A 中的 RUN₁ 和 RUN₂），或纯静息态记录（RUN₃）。会话片段指会话记录中的一个子部分，即受试者观看电影或处于静息状态的时间段。上述信息将在后续描述的 ISFC 计算中使用（见步骤 5.1.8），以限制不同会话片段类型带来的混杂影响。
   6. 在专用的可编辑文本窗口中输入与自助抽样（bootstrapping）相关的参数：执行 ISFC 计算的自助抽样折数，以及每次 ISFC 计算中构成参考组的受试者数量。
      注意：在我们的分析中，使用了 250 次自助抽样折数，参考组包含 6 名受试者。
   7. 在 时间参数 部分的专用可编辑文本窗口中，输入待分析时间序列的子时间段信息。需提供起始索引和终止索引（以 TR 为单位）。若要分析整个记录时长，请将起始索引设为 1，终止索引设为样本总数。
   8. 点击 绘图 按钮以执行 ISFC 计算。显示内容将随时间逐步更新，并同步显示已完成的自助抽样折数。对于区域对 (i,j) 和滑动窗索引 τ，ISFC 的计算方式为：在长度为 W 的滑动窗内，会话片段 s 与参考组中所有会话片段之间的交叉相关性的平均值；记参考组为 Ψ，其包含的受试者数为 N_Ψ，令 x_i^[s](t) 表示在时间 t 时会话片段 s 中区域 i 的时间序列，则 ISFC 估计值由下式给出：
      
      ISFC 测量在指定数量的自助抽样折数上进行计算，且在每一折中使用选定数量的会话片段作为参考组（见步骤 5.1.6）。若包含多个会话片段亚型，则参考组始终由混合的亚型样本组成。每个会话片段的最终输出是其未被包含在参考组中的所有折次中 ISFC 的平均值。
      注意：参考组是指在每次自助抽样折次中，会话片段 s 的功能时间序列所比较的对象集合。为了使结果对异常值更具鲁棒性，ISFC 在不同的参考组（即不同的会话片段子集）上多次计算。重要的是，采集时间 t 与滑动窗索引 τ 并不一致，因为后者基于 W 个数据点的集合计算，并依赖窗口步长大小获得连续估计值。该自助抽样过程借鉴了 Byrge 等人先前的一项研究³³。
   9. 若要保存输出结果以供后续步骤使用，请在专用可编辑文本窗口中输入保存名称，然后点击 保存 按钮。若要清空窗口内容，请点击 清除 按钮。

补充图3：首个与ISFC相关的图形用户界面窗口的示例截图。（上图）各被考虑会话中计算其ISFC测量值（即未被选为参考组成员）的频率示意图。（下图）在一名代表性受试者上，针对五十个示例连接计算得到的ISFC时间进程，这些连接是根据其在时间上绝对ISFC值总和最大而选出的。请点击此处查看该图的放大版本。

2. 在 MATLAB 终端中输入 JOVE_GUI4 以打开第二个与 ISFC 相关的图形用户界面窗口。
   1. 点击 加载 ISFC 数据 并选择在步骤 5.1 中生成的与刺激相关的 ISFC 输出文件
   2. 点击 加载空数据 并选择, 根据所使用的零假设数据生成方案，选择步骤5.1中生成的静息态ISFC，或相位随机化刺激相关ISFC输出文件
   3. 点击 加载编码本 并选择在步骤 4.3 中创建的码本文件。
   4. 在专用的可编辑文本窗口中输入数据的重复时间（TR，单位为秒）。
   5. 在专用的可编辑文本窗口中输入步骤5.1计算中使用的滑动窗口参数（窗口大小和步长，以TR为单位）。
   6. 在专用的可编辑文本窗口中输入（以百分比表示）用于对ISFC时间进程进行阈值化处理的α值，以突出显示显著性变化。
      注意： 此处及其他地方，当提及α值为2.5%时，表示当某一数值低于2.5%时即达到显著性水平^th 百分位，或大于97.5百分位^th 零假设数据的百分位数。在我们的分析中，共有 5,762 个静息态数据点可供使用，并选择 α 值为 10^-4这意味着我们希望有0.01%的数据样本大于或等于所选阈值，超过该阈值的ISFC偏移将被视为具有显著性。作为对比，经Bonferroni校正后所需的α水平为0.05/44,551 = 1.12 × 10⁻⁶^-6，而根据我们的数据量（n 个样本）所能达到的最严格的 α 水平为 .
   7. 点击 绘图 点击按钮执行ISFC阈值化过程，该过程会针对给定连接聚合所有可用的空ISFC测量值，以构建空分布，随后根据所选α值对与刺激相关的ISFC测量值进行阈值化处理。当与刺激相关的ISFC值在统计学上显著超过空分布时，对应的时间点将被标记为-1（表示ISFC显著降低）或+1（表示ISFC显著升高）。
      注意： 阈值化过程借鉴了 Betzel 等人关于静息态动态功能连接的研究工作。²³.
   8. 要可视化不同时间点的ISFC空间模式，请拖动ISFC轨迹图下方的滑块。

补充图 4：第二个与 ISFC 相关的图形用户界面窗口的示例截图。（左上图）在代表性受试者上，计算了三个示例连接的 ISFC 时间序列，这些连接被选为表现出最显著 ISFC 偏离的连接，并显示其相应的计算显著性阈值（水平线）。（左下图）对于相同的连接，显示了跨受试者平均的偏离时间序列，误差量以双尾 95% 置信区间表示。（右图）在 ISFC 图和偏离图中由一条黑色垂直线指示的选定时间点上的空间 ISFC 模式（跨受试者平均的 ISFC 偏离）。正向 ISFC 偏离以黄色表示，负向偏离以粉色表示。节点的大小和颜色编码与其度值成正比。请点击此处查看该图的放大版本。

结果

本研究纳入了15名正常发育（TD）的受试者，均已获得书面知情同意。所有受试者均为右利手男性（年龄为23.42 ± 7.8岁）。所采用的实验范式是一部面向青少年的视听科学纪录片，内容涉及日光暴露的危害。该纪录片包含丰富的视觉、听觉和社会性刺激，可通过以下网址观看：https://miplab.epfl.ch/index.php/miplife/research/supplement-asd-study。

每位受试者采集了两个会话数据（RUN₁ 和 RUN₂），其中播放评估影片的时间为 5 至 353 秒（持续 5.8 分钟）。随后还采集了一段静息态数据，时间为 386 至 678 秒（持续 4.9 分钟）。此外，每位受试者（除一名因幽闭恐惧症未完成者外）还单独采集了一次仅包含静息态的会话（RUN₃），持续 310 秒（5.2 分钟）。示例电影场景及数据采集时间点总结于图 1A中。需要强调的是，该采集方案并非最优，因为电影刺激后立即进行的静息态记录可能部分受到溢出效应的影响²⁷；本研究利用这些数据以获得足够数量的样本用于统计阈值分析，但应尽可能避免此类设计。

我们排除了所有在0.5 mm阈值下被剔除帧数超过10%的扫描会话，并采用Craddock 等人³⁴（双层时间相关算法）的脑区分割方法，生成区域时间序列，共涵盖299个不同的脑区。

ISFC 分别在以下三部分数据上独立计算：(1) RUN1 和 RUN₂ 中观看电影的子阶段，(2) RUN₁ 和 RUN₂ 中静息态的子阶段，以及 (3) 静息态 RUN₃ 的记录数据。主要结果采用时间窗长度 W = 10 TR，并与较短的 W = 5 TR 进行比较。步长始终设为 1 TR。重采样（Bootstrap）分析共进行 250 次迭代，每次参考组包含 6 个扫描会话片段。

图1B展示了在 W = 10 TR 和 W = 5 TR 条件下生成的ISFC时间过程，涉及三个不同且具代表性的功能连接：连接1包含左侧下顶叶区域（与运动物体预期相关，MNI坐标：41,9,32）³⁵和右侧额盖区域（与反应抑制相关，-34,-52,45）³⁶。后一脑区还分别参与了连接2和连接3，其中连接2涉及与感觉协调相关的脑区（54,6,34）³⁷，连接3涉及与词语意义加工相关的脑区（6,62,9）³⁸。

不同时间窗长度的比较显示，在 W = 5 TR 条件下，无论是观看电影还是静息态记录，被试间的时间方差整体上均大于 W = 10 TR 条件，这是滑动时间窗分析中已知的现象³⁹。对于连接1，无论时间窗长度如何，观看电影记录中的一个局部时间段（约55秒处）均表现出被试间强烈的、同步的ISFC升高，其值远超静息态情况下的数值范围。因此，我们期望通过阈值化方法将该时间段识别为一个显著的ISFC瞬态事件。

对于连接2，我们观察到类似的时间动态，但对于W = 5 TR的情况，由于滑动窗口方法引入的噪声较大，其信号变化较难与静息态时间过程区分开来。至于连接3，它反映了一种对影片刺激无明显响应的情况，因此观影状态与静息态时间过程的波动模式相似。在本分析阶段的预期结果是：部分连接表现出明显的刺激诱发重组，而另一些连接则无响应，两者共同构成混合结果。

图1：数据采集时间安排及ISFC时间序列示例。（A）受试者观看的电影包含多种社交情境（示例图像1和4）、配有彩色图示的科学解释内容（示例图像2和5）以及自然景观场景（示例图像3）。每位受试者采集了三个扫描会话：其中两个会话（RUN₁ 和 RUN₂）包含电影刺激阶段（从5至353秒，绿色高亮）及随后的静息态阶段（从386至678秒，黄色显示）；另一个会话（RUN₃）仅包含静息态记录（持续310秒，橙色显示）。（B）在三个代表性功能连接（C1、C2和C3，分别对应深绿/红色、浅绿/橙色及青绿/黄色曲线）上，电影观看期间（冷色调）与静息态期间（暖色调）ISFC随时间的变化情况。当时间窗W = 10个TR时（左图），电影观看期间ISFC的变化比W = 5个TR时（右图）更为显著。每条曲线代表一个会话中的ISFC时间序列。本图部分改编自Bolton 等人²⁵。请点击此处查看此图的放大版本。

图2A展示了对ISFC时间序列进行统计阈值处理后的结果，涉及与上述相同的三个连接。时间序列值为1表示所有受试者在同一时间点均出现相同的ISFC升高；值为0表示没有任何受试者出现显著的ISFC变化；值为-1表示所有受试者同步出现ISFC降低。与之前相同，我们对比了W = 5 TR与W = 10 TR的情况，同时突出了两种α值情形：α = 0.01% 和 α = 5%。

与上述观察结果一致，较短的时间窗长度会减少所提取的显著ISFC变化的数量。对于连接1，W = 5 TR和W = 10 TR均识别出同一特定时刻（t = 55 s）出现明显的ISFC增强。考虑到约5秒的血流动力学延迟，这一时刻对应于电影中彩色线条向一个玩偶延伸并突然在其前方停止的片段（46–49 s），这与相关脑区在运动物体预期和反应抑制中的作用相符^35,36。

当将α从0.01%增加到5%时，可以观察到检测到的ISFC瞬态特异性显著降低，可能包含大量假阳性结果，并且预期表现出更少的时间同步性。

从另一个角度观察数据，图2B展示了在 t = 55 s 时显著的ISFC变化的全脑空间图谱。可以看出，对电影场景的响应远远超出了此处描述的示例连接范围。

图2：ISFC模式的时间与空间快照 (A) 三类代表性连接（C1、C2 和 C3，分别对应深绿色、浅绿色和青绿色曲线）在受试者间的 ISFC 瞬态时间进程的平均结果。浅灰色区域标示出引发 ISFC 变化的电影片段，并以示例图像展示。当 W = 10 TR 时（左图列），ISFC 变化的检测强度高于 W = 5 TR 时（右图列）。当 α = 0.01% 时（上图行），对局部电影线索的特异性高于 α = 5% 时（下图行）。每条曲线代表一次实验会话的 ISFC 瞬态时间进程，双尾 95% 置信区间作为误差度量显示。B当 W = 10 TR 且 α = 0.01% 时，在 t = 55 s 时（即 C1 的 ISFC 瞬态峰值时刻），ISFC 瞬态呈现出清晰且局限的空间分布模式；而当 W = 5 TR 且 α = 5% 时，此时发生显著 ISFC 变化的连接数量则多得多。需注意，我们假设所述时间关系中存在约 5 s 的血流动力学延迟（即, 此处的55秒值对应于50秒时的视频刺激。该图部分修改自 Bolton 等人。²⁵. 请点击此处查看此图的放大版本。

讨论

在一组健康受试者的数据集中，我们展示了功能连接（FC）的同步跨被试增强与减弱（即ISFC瞬态）如何在时间上与局部化的电影线索相匹配，从而提供超越静态描述的信息。尽管使用跨被试相关性度量能够使分析聚焦于刺激驱动的功能重组，但也必须意识到这种方法将研究结果局限于研究人群中普遍存在的效应：因此，与前额叶加工相比，低级感觉加工的表征可能被过度放大⁴⁰。为克服这一局限，目前正发展一些新方法，这些方法还能够提取在不同被试间变异最显著的脑区⁴¹。

所提出方法的另一个局限性源于滑动窗口的设计，与逐帧分析方法相比，ISFC瞬态时间过程的时间分辨率有所降低¹⁵。正如我们所展示的，需要在窗口长度足够短以充分解析动态ISFC重组，以及窗口足够大以获得稳健估计之间进行权衡。在我们的框架中，有两个关键步骤可确保提取出的ISFC瞬态真实反映了连接性发生的实际变化：第一，使用窗口长度的倒数对区域时间序列进行高通滤波³²；第二，利用静息态ISFC数据生成相关的零分布，并且其采集参数与刺激相关数据完全相同。当然，后者也要求更长的总体采集时间，以便在刺激相关扫描之外额外获取静息态数据。作为一种避免额外静息态记录的替代方法，我们还提供了直接从刺激相关时间序列生成相位随机化数据的可能性，这种方法在动态功能连接分析中常被采用^23,24。对部分扫描会话的进一步评估表明，尽管静息态零分布方法更为保守，因而更不易产生假阳性结果，但两种方案在ISFC偏离事件检测的整体模式上是相似的（见补充图5）。

补充图5：不同零模型数据生成方法下ISFC偏移的检测。针对静息态（左列，蓝色图）或相位随机化（右列，红色图）的零模型数据生成方法，提取出的跨连接ISFC偏移百分比。底部图为源自前三个脑区的连接的放大图。误差条表示受试者间的标准差。请点击此处查看该图的放大版本。

静息态数据采集的时长实际上与分析中的一个关键参数——α值相关。如上所述，若选择过于宽松的α值，将导致检测到的ISFC瞬变信号中出现大量假阳性结果。可用的静息态数据量越大，可实现的假阳性率就越严格，因为阈值的设定可以依据零分布中更为极端的数值。作为参考，在本研究中采用n = 299个图谱区域，并拥有共计5,762个静息态数据点的情况下，我们能够达到的最佳α值接近0.01%（数学细节参见步骤5.2.6）。

与任何fMRI分析相关的另一个关键点在于，需严格去除分析数据中可能存在的与运动相关的伪影^30,42。特别是，如果希望将本研究所介绍的流程应用于在扫描过程中表现出明显运动的疾病人群，我们建议除了在统计分析中将运动变量作为协变量纳入之外，还应增加额外的预处理步骤，例如小波去噪⁴³或ICA-AROMA⁴⁴。组间比较（例如，比较健康组与疾病组之间的ISFC瞬态变化）可通过对两组目标人群并行执行所述方法来实现（参见Bolton等²⁵关于自闭症谱系障碍诊断人群的示例）。然而，组间差异可能出现在两种不同情形下：（1）某一组中缺乏ISFC变化，或（2）该组中变化模式更具异质性。为了区分这两种因素，应对疾病组再次运行该流程，在自举（bootstrapping）步骤中使用健康受试者组作为参考组。前一种情况仍将表现为无响应，而后一种情况则不会。

除了我们在此处描述的内容之外，所介绍的方法还开辟了颇具前景的未来研究方向：从分析角度来看，ISFC瞬时图可被视为脑网络图，从中可提取量化脑功能连接的指标⁴⁵，或可通过聚类方法提取动态ISFC状态，并对其空间和时间特征进行评估^17,46。此外，还可以考虑使用比皮尔逊相关系数更复杂的连接性测量工具，以揭示功能连接（FC）更细微的特征^47,48。

从实验角度来看，将我们的分析流程应用于更广泛的范式是一个具有前景的方向：例如，除了本文所研究的电影之外，还可以考虑使用一段音乐⁴⁹或一个叙事性故事^13,50作为时间锁定的刺激材料。或者，甚至可以通过超扫描技术⁵¹来探究自然情境下的社会交流过程^52,53。

材料

本文使用的材料清单

姓名	公司	目录编号	评论
Freesurfer 版本 6.0	计算神经影像学实验室，马萨诸塞州波士顿马丁诺斯生物医学影像中心，美国	https://surfer.nmr.mgh.harvard.edu/fswiki/DownloadAndInstall	一个与 MATLAB 兼容的工具箱，可用于对功能磁共振成像数据执行各种处理、可视化和分析步骤
MATLAB_R2017a	MathWorks	https://ch.mathworks.com/downloads/	可运行的 MATLAB 计算软件版本（应使用 2014a 或更新版本）
统计参数映射版本 12.0（SPM12）	惠康信托神经影像中心，伦敦大学学院，英国伦敦	https://www.fil.ion.ucl.ac.uk/spm/software/spm12/	一个与 MATLAB 兼容的工具箱，可用于对功能磁共振成像数据进行统计分析
Tim-Trio 3 T MRI 扫描仪	西门子	https://www.healthcare.siemens.ch/magnetic-resonance-imaging/for-installed-base-business-only-do-not-publish/magnetom-trio-tim	用于记录受试者脑功能活动的磁共振成像扫描仪（场强为 3 T）