静息态fNIRS中的全局生理信号校正

导读

意义：功能性近红外光谱（fNIRS）数据的静息态功能连接（RSFC）分析揭示了整个大脑的皮层连接和网络。fNIRS信号中的运动伪影和系统生理学提出了独特的分析挑战，目前已有研究探索了控制系统生理噪声的方法。但是这些相同的方法在应用于静息态fNIRS（RS-fNIRS）数据时是否需要修改仍不清楚。

目的：本研究系统地考察了几种RSFC分析管道的敏感性和特异性，以确定校正RS-fNIRS数据中全局系统生理信号的最佳方法。

方法：使用数值模拟的RS-fNIRS数据，比较了几种连通性分析管道的真阳性率和假阳性率。使用受试者工作特征分析对其性能进行评分。管道包括偏相关和多变量格兰杰因果关系，有/没有短距测量，以及包括非传统零滞后交叉项的修正多变量因果关系模型。研究者还考察了预白化和稳健统计估计对性能的影响。

结果：与以往关于双变量相关模型的研究一致，本研究结果表明稳健统计和预白化是校正fNIRS时间序列中运动伪影和自相关的有效方法。此外，当两个感兴趣的信号同步波动时，使用从短距fNIRS通道提取的主成分作为偏相关模型的一部分进行预滤波，在减少由于共有的系统生理学引起的伪相关方面最有效。然而，当信号之间存在时间滞后时，结合短距通道的多变量格兰杰因果检验效果更好。由于不知道实验数据中是否存在这样的滞后，研究者提出了一种修正版本的格兰杰因果关系，其中包括非传统的零滞后项作为折中解决方案。

结论：建议结合预白化、稳健统计方法和处理管道中的偏相关性，以减少自相关、运动伪影和全局生理噪声，从而获得具有统计有效性的RS-fNIRS连通性指标。未来的研究应使用人类数据来验证这些方法的有效性。

前言

自发或静息态功能连接（RSFC）已成为表征人脑网络结构的重要工具。Biswal等人首先指出，功能性磁共振成像（fMRI）血氧水平依赖性（BOLD）信号中的同步自发波动被组织成不同的网络，这些网络在被试之间保持不变。25年后，这些结果在数百项研究中被复制，并扩展到人类大脑功能、发育和疾病的各个方面。静息态网络也使用其他神经成像方法进行了观察，包括EEG和MEG，以及fNIRS。静息态范式的一个主要优点是不需要执行任务，这在婴儿、老年人、患者等人群的研究中特别有用。

目前已经开发了几种分析方法来量化静息态连通性。对于fMRI和fNIRS，这些方法通常关注大脑区域之间的零滞后关系。零滞后关系是指两个信号同时发生信号波动（在测量采样时间内）的关系。估计功能连通性的最常用方法可能是计算两个脑区时间序列之间的Pearson相关系数。有几种方法可以在时域中表征RSFC，包括（i）基于种子的连接，（ii）成对连接，以及（iii）独立成分分析（ICA）。在基于种子的连接中，选择一个种子区域，并评估其与所有其他脑区域的连通性。在成对连接中，计算所有可能通道对之间的相关性。ICA将fNIRS通道数据的协方差分解为几个统计上独立的时空成分，其中一些成分表示每个被试的大脑功能网络，或者表示跨被试连接的功能网络。

格兰杰因果关系是“有效连接”方法的一个例子，该方法依赖于fNIRS通道在时域上的滞后关系。另一种常用于模拟fNIRS中功能连通性的方法是小波变换相干（WTC），该方法在时频域中模拟通道之间的关系。WTC依赖于连续小波变换（CWT），CWT将时间序列数据转换为时频域。WTC估计为时间序列在时频域中的互相关。

功能性近红外光谱

fNIRS是一种光学成像方法，使用近红外波长范围（650至950nm）的光来估计氧合血红蛋白（HbO）和脱氧血红蛋白（HbR）浓度的时间变化。fNIRS测量大脑的血流动力学变化，反映血流量、体积和血氧合的变化，因此，它是一种间接测量神经活动的方法，就像fMRI BOLD信号一样。与fMRI相比，fNIRS具有一些优势，包括成本更低，对头部运动的弹性更大，以及更高的数据采样率。fNIRS也可以便携，能够执行更灵活或生态有效性的任务。然而，fNIRS的缺点是其空间分辨率较低，深入大脑的深度有限，相对于大脑本身，对头部表层的敏感性增加。对于fNIRS测量来说，后一个问题尤其令人困惑，因为头皮的系统生理信号，包括心脏、呼吸和血压波动，通常会导致空间上的全局噪声，其振幅通常大于感兴趣的潜在大脑信号。因此，在使用fNIRS的RSFC实验中，将表层生理噪声和真正的潜在神经相关脑网络进行分离是具有挑战性的。

Obrig等人首次证明了fNIRS在捕获极低频和低频振荡方面的可行性，这些振荡有助于RSFC的研究。之后，几项用于考察RSFC可行性的早期研究使用fNIRS和扩散光学层析成像（DOT）发现了可靠的和预期的视觉、感觉运动和语言网络。从那时起，研究人员使用静息态fNIRS（RS-fNIRS）和DOT研究婴儿发育、前额叶皮层的性别差异、早产、语言网络、自闭症谱系障碍、情感障碍和衰老。检查RSFC的重测信度发现，在从一小时到一周的多个sessions中，个体和群体水平的RSFC图之间存在良好到极好的组内相关性。

RS-fNIRS信号的统计特性

鉴于RS-fNIRS的信号和噪声特性，与RS-fMRI相比，RS-fNIRS具有一些独特的挑战。由于fNIRS信号的采集方式，它对头皮的非神经脑外信号变化很敏感。RS-fNIRS信号由三部分组成：（i）低频神经振荡（感兴趣信号），（ii）源自大脑的脑内生理噪声，以及（iii）源自头皮和非脑组织的脑外生理噪声。鉴于静息态范式的性质，所有信号都是非诱发的。因此，fNIRS信号同时受到脑外信号变化和非神经脑信号变化的影响，包括心脏（约1-1.2Hz）、呼吸（0.3-0.6Hz）和血压/Mayer波（0.1Hz）波动引起的系统生理噪声。生理噪声对RS-fNIRS信号有两种影响：①诱导时间自相关和②增加大脑通道之间的空间协方差。

先前的研究表明，高采样率和系统生理信号(心率、呼吸和血压)，结合潜在的血流动力学反应，可以共同导致时序相关，从而在RS-fNIRS时间序列中会产生自相关。去除时间序列中的自相关以获得统计有效的关系估计，可以使用预白化的方法。预白化是通过去除时间序列中的序列相关来实现的。自回归(AR)模型可用于数据预白化和去除序列相关。具体来说，对于fNIRS，使用预白化需要使用p阶AR模型来模拟当前时间点作为前一个时间点的函数。AR模型中去除时间相关性后的残差项被用来计算通道之间的Pearson相关系数，而不是使用自相关通道数据。当将预白化纳入到预处理管道中用于功能连接估计时，可以减少I型错误。尽管对fNIRS的时间自相关进行了校正，但如果没有考虑由于系统生理机能导致的空间协方差增加，这也可能会导致FDR的增加。

由于系统生理信号被认为在整个大脑中是相对均匀的，所以由所有通道共有的系统波动引起的信号变化可以增加通道中的空间协方差。图1显示，不校正fNIRS数据中增加的空间协方差可能会导致I型误差增加和特异性降低，因为基于地面真值的连接路径远比预期的多。

图1.连接图显示，在模拟的fNIRS数据中，加入模拟AR平滑和共享协方差的生理噪声，Pearson相关系数、AR相关和AR偏相关可以有效地恢复真实的连接。

方法

为了表征每种分析方法的敏感性和特异性，本研究根据具有已知地面真值的数值模拟数据，探索了模型在几种模拟情况中的性能，从简单正态分布随机噪声（其中本研究模型的所有统计假设都有效）到具有时间自相关噪声的生理信号，再到全局空间生理噪声。最后，加入运动伪影。每个数据模拟都代表了一个更具挑战性的统计问题，通过模拟数据中的真实和零相关性，使用受试者工作特征（ROC）分析来检验这个问题。所有方法均在MATLAB（R2019b）的NIRS AnalyzIR工具箱中实现。

数据模拟

采用具有指定空间协方差矩阵的多变量正态分布方法，对由16个长距通道（27mm）和9个短距通道（7.1mm）组成的探测器进行了fNIRS数据模拟。模拟了近红外波长为690和830nm的数据。这种已知的空间协方差在数据中产生潜在的“真实”相关性，在此基础上评估各种算法。这些真实相关性是通过通道之间的零滞后（无时间偏移）或时间滞后（数学因果模型）产生的。仅针对长距通道在通道空间中模拟目标真实相关性。为了给模型添加额外的噪声，使用半无限均匀平板模型来计算光学正向模型。使用光学正向模型将这些体素投影到通道空间，并添加到基本真值相关信号。

对于每个模拟，以4Hz采样率生成300 s的静息态数据。为了生成ROC曲线，每个测试总共生成了200个数据集。在ROC分析中，在探测器中任何两个通道之间的所有可能连接中，模拟了10%左右的可能连接路径的真实连接，确保全局系统噪声是主要的空间噪声特征。因此，用于生成模拟数据的协方差矩阵中只有10%的非对角元素为非零。由于时间自相关、全局系统生理信号或运动伪影，其余90%路径之间的任何显著连接都是假阳性。为了生成真和假阳性估计，从总邻接矩阵中随机抽取相等数量的真和空连接。

fNIRS实验数据

采集

研究者还测试了本研究的功能连接方法对实验静息态数据的有效性和适用性。使用商业NIRScout-2（NIRx Medizintechnik GmbH，柏林，德国）连续波fNIRS系统记录NIRS数据。共有50个通道（42个长距通道，8个短距通道）分布在双侧额叶和感觉运动皮层（见图2）。如图2所示，长距通道由放置在头皮上的16个源光电二极管(橙色圆圈)和13个探测器(蓝色矩形)组成。8个探测器（绿色菱形）放置在头皮的8个短距通道上。蓝色实线表示长距通道，绿色虚线表示短距通道。以7.8125Hz的采样率记录两个波长（760nm和850nm）的数据。定位头帽后，使用NIRx Aurora软件优化信号质量。

图2.fNIRS数据采集示意图。

参与者和任务

约24名健康被试参与了实验（9名男性和15名女性；平均年龄27.4岁，SD=8.7岁；23名右利手）。告知被试实验情况，并获得书面同意。这项研究得到了匹兹堡大学机构审查委员会的批准。每个被试进行一次静息态扫描和其他任务扫描。指示被试在静息态扫描期间尽量减少身体运动，保持坐姿放松5分钟，而不使用任何脑力劳动。

数据预处理

将模拟的fNIRS信号强度转换为光密度的变化，然后使用修正的Beer-Lambert定律（MBLL）将光密度单位转换为HbO和HbR浓度，两种波长的微分路径长度因子为6，部分体积校正均为60。除了对实验fNIRS时间序列进行降采样外，以与模拟数据类似的方式对数据进行预处理。使用MBLL将原始信号强度转换为HbO和HbR浓度后，将实验静息态fNIRS数据从7.8125Hz降采样至4Hz，以匹配模拟数据的采样率。在本文的其余部分，分析的重点在于HbO浓度的变化，结果也应适用于HbR。

图3显示了三种情况下模拟数据的fNIRS时间序列示例，（i）仅时间自相关，（ii）具有时间自相关和全局系统生理信号，以及（iii）具有时间自相关、全局系统生理和运动伪影。图3(左)显示了在使用MBLL将原始模拟fNIRS信号转换为HbO信号变化后，三种情况下的两个通道的时间轨迹。图3(右)显示了通道的相关结构。正如预期的那样，随着全局系统生理信号的引入，这些通道中的协方差增加。随着在模拟数据中引入运动伪影，随着通道之间的相关值接近1，协方差进一步急剧增加。

图3.显示了模拟数据中两个代表性通道HbO变化的时间轨迹以及通道的归一化协方差矩阵。

预白化

预白化通常使用AR模型来去除时间序列中的自相关，并对信号的频率内容进行白化。预白化大大降低了连接模型中的FDR。这里，使用AR模型实现预白化，如下所示：

fNIRS时间进程可表示为：

其中εt∈N(0，σ2)，yt是fNIRS的时间序列，ai是AR参数，p是模型阶数，由贝叶斯信息准则（BIC）确定，εt是新息项，用于计算连通性指标，而不是时间序列本身。图4显示了三组模拟数据预白化前后代表性通道的（a）功率谱密度和（b）自相关函数图，包括：（i）时间自相关，（ii）时间自相关和全局系统生理信号，（iii）时间自相关、全局系统生理和运动伪影。正如预期的那样，预白化步骤使功率谱密度图变得平坦，并降低了fNIRS时间序列中的时间自相关。

图4.显示了三组模拟数据预白化前后代表性通道的(a)功率谱密度和(b)自相关函数图。

稳健方法

在对时间序列进行预白化后，新息项包含在该时间点添加的信息。由于运动伪影被视为异常值，它们可能会影响正常性假设。通常，这些高影响点会对连通性指标的估计产生重大影响。如前所述，运动伪影对连接性指标的影响不同，这取决于运动是跨多个通道共享还是仅限于几个通道。跨多个通道共享的运动损坏时间点可以通过增加它们之间的协方差来影响功能连接。因此，为了识别统计异常值，研究者使用了之前在Santosa等人研究中引入的联合加权函数的多变量版本。联合加权函数可以通过计算时间序列的几何长度来识别异常值，并根据它们偏离平均值的方式对它们进行适当的降权，将运动损坏的异常值分配较低的权重，使它们对相关和回归方法的影响较低。

根据所有通道（包括短距通道）的新息项计算几何长度函数rt，如下所示：

其中εi，t表示时间点t的第i个通道的新息项，c是通道总数。

根据离群值偏离平均值的程度对离群值进行加权的加权函数(S)，可使用Tukey双平方函数的平方根，对所有时间点进行如下计算：

这里，σ是分布的标准偏差，根据中值绝对偏差（MAD）和正态分布标准偏差（σ）之间的关系估计(σ≈ 1.4826*MAD）。使用MAD是因为它是更稳健的分散统计。常数κ的值为4.685，这对于没有任何异常值的标准正态分布提供了95%的效率。具有对角元素的加权矩阵S表示所有时间点的加权函数可与所有时间序列的新息项相乘，从而降低正态分布统计异常点的权值。

其中，Yc是所有通道重加权后的新息项矩阵，Ec是预白化后获得的新息项矩阵。

因此，为了减少运动伪影的影响，对于回归模型，不使用线性模型和普通最小二乘估计参数，而是使用稳健统计，使用加权最小二乘（WLS）来降低从学生化残差计算出的异常值的权重。权重是使用Tukey的双平方函数计算的。然后将计算出的权重应用于原始数据，并重新估计残差。重复该过程直到收敛。稳健回归在内置的MATLAB函数“robustfit”中实现。MATLAB函数robustfit使用50次迭代或回归系数值无可检测数值变化的收敛标准来评估系数的收敛性。由于多变量格兰杰因果关系（MVGC）和修正的MVGC在50次迭代中不收敛，研究者将最大迭代限制更改为150次迭代，以确保回归参数收敛。

使用WLS回归估计的参数可用如下公式获得：

这里，是参数估计，S是加权矩阵。

用于减少全局信号变化的基于PCA的回归

如果探针广泛分布在大脑中，那么跨多个通道共享的大信号方差可能主要反映生理噪声。因此，PCA的前几个成分可用于识别由于生理噪声引起的全局信号变化，并从通道数据中回归其影响。数据协方差矩阵可以使用奇异值分解法进行分解

其中C是给定C=YT·Y的数据协方差矩阵，U是包含左奇异向量的矩阵。∑是包含奇异值的对角矩阵，V是包含右奇异向量的矩阵。通常，使用解释最大方差的前几个成分。研究者使用最小数量的主成分来解释至少90%的通道协方差。对长距通道、短距通道或两者进行协方差分解可以提取全局信号。这三种选项都在测试中进行了检查。经PCA回归后的通道数据如下：

这里，y是新息项或通道的原始数据，取决于分析通道。Z是包含前几个成分作为列的矩阵。

功能连接指标

本研究考察了几种估计功能连接的方法。在这里，重点关注时域中的方法，包括相关性、偏相关和格兰杰因果关系的几种变式。时频方法，如小波相干，由于两个原因没有在本研究中进行检查。首先，在研究者之前的工作中发现，在双变量连接性分析中，小波相干优于ROC分析中的标准Pearson相关，但性能不如预白化的稳健相关模型。除非标准小波相干模型被修改为包括高阶AR项和稳健（异常值拒绝）方法，否则小波相干也具有较高的假阳性率和不受控的I型误差。经过这些修改发现，与预白化稳健模型相比，小波相干没有任何优势。第二，小波相干的多变量或偏相关版本尚未描述，需要开发一个新的统计模型。

基于相关性的分析模型

Pearson相关系数。两个时间序列x和y之间的Pearson相关系数计算如下：

其中σx，y是时间序列x和y的协方差，σx和σy分别是x和y标准差。标准Pearson相关性并不是一个稳健的统计估计，Santosa等人在fNIRS的背景下提出并描述了若干稳健变量。稳健相关性计算为回归模型y对x和x对y的标准回归系数的几何平均值，以减少x或y中的异常值，但这仍然会影响相关估计。步骤如下：

基于上述回归模型中的学生化残差估计加权矩阵，加权矩阵Sx→y在上述回归模型的两侧相乘，并且更新参数估计。

同样，可以在y上拟合x，

加权矩阵Sy→x根据学生化残差估计，以降低x中的异常值。Sy→x在上述回归模型的两侧相乘，并更新参数估计。

两个回归模型交替求解，直到参数βy→x和βx→y收敛。使用以下公式估计最终稳健相关系数：

稳健相关系数估计的符号由回归系数的符号确定。有效自由度为

。通过降低数据点的权重，特别是对异常点的权重，降低了有效自由度。

对于标准和稳健相关估计，在给定零假设（p值）的情况下，观察数据的概率可以通过使用T统计量估计

其中df为标准Pearson相关的有效自由度，值为n−2，r为Pearson相关系数。这里n是时间点的数量。

自回归相关。AR相关性是基于预白化时间序列计算的Pearson相关系数。如Santosa等人所述，为每个数据通道估计p阶AR模型，并用于预白信号。然后在相关分析中使用AR模型中的新息项（例如，在每个时间步长添加的新信息）。研究者使用的最大模型阶数为40（10×采样率），高于模拟fNIRS数据时使用的AR滤波器。根据经验观察，在实验静息态fNIRS数据中，模型阶数通常＜10·Fs（此处Fs为采样频率）。

AR偏相关。偏相关涉及一个两步过程，其中首先将其他通道从感兴趣的数据中投影出来，然后对残差执行相关。使用基于PCA回归的投影来部分消除其他通道的影响。PCA回归用于避免模型的共线不稳定性。使用最小数量的特征向量来解释至少90%的空间协方差。对于两个通道i和j，如果预白化后的新息项由εi和εj表示，则白化PCA滤波时间序列计算为

现在，当从通道中去除全局信号时，AR偏相关可以计算为两个残差之间的Pearson相关性。基于保留的主成分数量进一步减少自由度。该模型的输入数据是预白化新息项时间序列，而不是通道数据，因此该方法被称为AR偏相关。

基于格兰杰因果关系的分析模型

多变量格兰杰因果关系。在控制其他时间序列的影响后，MVGC检验两个时间序列之间的滞后关系。这本质上是两个回归模型的模型拟合检验。可以使用多变量AR模型（MVAR）控制通道Z后检验x“格兰杰原因”是否为y。对于通道时间序列，可以将当前时间点yt建模为限制性模型中y和其他通道Z之前滞后的函数：

（限制性模型）

在非限制性模型中，添加x的滞后项，以模拟y中的当前时间点：

（非限制性模型）

这里，p是模型中包含的滞后项的数量，使用BIC的拟合测量确定。为了检验因果关系，进行了以下统计检验，与没有x（限制性模型）的模型相比，在x（非限制性模型）中添加额外滞后项是否提高了yt的可预测性。如果x的过去包含关于当前时间点yt的预测信息，则即使在控制丢失的自由度之后，yt预测中的误差也会得到改善。这在统计上可以评估为：

零假设：x不是y的格兰杰原因

H0: β1=β2=β3=……βp=0

备择假设：x是y的格兰杰原因。

HA：至少有一个βi为非零

因为有嵌套模型，所以可以使用嵌套F统计量来检验格兰杰因果关系

这里，SSER和SSEU分别是限制性和非限制性模型的平方误差之和，n是时间点的数量。因此，可以在自由度p和[n−(2p+1)]处计算F统计量的p值。格兰杰因果指标G也可以用于评估滞后关系的强度，如下所示：

在双变量格兰杰因果关系上使用MVGC的目的是控制包括系统生理信号在内的第三变量的影响。在MVGC模型中添加多个通道的滞后项，容易使回归模型无法求解，因为变量的数量将超过时间点的数量。因此，与使用通道数据本身不同，首先执行PCA分析，并像之前使用偏相关方法一样，包含解释至少90%方差的最少成分。然后，在限制性和非限制性模型中包括成分和时滞成分，而不是通道数据。

修正的多变量格兰杰因果关系。传统的MVGC仅包括方程右侧的y、x和Z的时间滞后历史。因此，MVGC仅对非零滞后关系建模。重要的是要认识到零滞后与采样率有关，在模拟fNIRS数据的4Hz采样时，一个采样（250ms）内发生的任何同步信号变化在数学上都是零滞后。虽然大脑中的神经信号传递相当快（100ms），但我们不希望假设所有的相关性都是零滞后，因而本研究引入了一个改进的MVGC模型，以建模通道之间的滞后和零滞后关系。在修正的MVGC中，将Z的零滞后项包括在限制性和非限制性模型中，以及通道x的零滞后项包括在非限制性模型中。因此，限制性模型可以写成：

（限制性模型）

类似地，非限制性模型可以写成：

（非限制性模型）

上述方程还包括Z和x的零滞后项。现在，为了检验格兰杰因果关系的统计显著性，检验了如下零假设和备择假设：

H0: β1=β2=β3=……βp=0

HA：至少有一个βi为非零

统计显著性采用嵌套F统计进行评估，格兰杰因果关系的计算与传统MVGC类似。该修正模型的一个警告是，假设不检验y和x之间的连接是否是由于零或非零滞后项。我们可以通过修改非限制性模型以包含x的非零滞后项来进一步检验这一点，使得非限制性模型和限制性模型之间的唯一区别是在x中包含零滞后项。

稳健的多变量格兰杰因果关系。最后，研究者开发了MVGC和修正的MVGC模型的稳健统计版本。该多变量模型中离群值的定义具有挑战性，因为如双变量相关模型，离群值可能存在于y、x和/或Z的时间过程中。在双变量相关中，使用离群值权重的并集（例如，y或x中的离群值）。然而结果发现，扩展到多个通道的数据通常会导致太多的数据被向下加权，因为任何一个通道中的异常点会导致所有通道的时间点被向下加权。特别是，MVGC中的时滞要求包含该列矩阵的所有滞后列也被向下加权。

为了降低上述回归模型中使用多个项的稳健回归方法的复杂性，研究者使用新息项，而不是单独使用通道数据。稳健MVGC的过程如下所述。首先，对通道数据（短距或长距）进行预白化，并计算新息项。使用从新息项计算出的加权矩阵对数据进行预加权。经过预白化和预加权后，其他通道（Z）的新息项包括在限制性模型的回归中，如下所示：

（限制性模型）

这里，εty是在预白化和预加权后的时间点t处通道y的新息项。εjz,t−i是时间点t-i处第j个成分的主成分权值。nPC表示使用的主成分的数量。

同样，新息项也用于非限制性模型，而不是通道数据本身：

（非限制性模型）。

代码实现

所有方法均在MATLAB中实现。利用MATLAB内置函数和NIRS Brain AnalyzIR工具箱对相关方法进行预处理和实现。使用自定义的MATLAB脚本实现了偏相关、格兰杰因果关系和稳健格兰杰因果等方法。本文中使用的工具箱和代码可在线访问：https://github.com/pradlanka/rsfc-fnirs获得。

结果

处理方法的比较

如预期那样，随着模拟的进展和噪声的增加，标准Pearson相关等基本模型很快就难以为继了。例如，Pearson相关假设独立的测量，这与加入的时间自相关噪声相违背。如图5(b)所示，与预白化和MVGC模型相比，不包括AR（预白化）滤波器的Pearson和稳健Pearson相关模型在ROC曲线上的性能显著更低。几乎所有其他方法的AUC接近1。此外，如图5(a)所示，Pearson和稳健Pearson模型具有较大的假阳性率。其他分析模型的I型误差图更接近理想（对角线）。AR相关和稳健AR相关都具有最接近理想的I型误差控制。AR偏相关模型略微高估了测试统计的显著性，而稳健的AR偏相关过校正导致了对测试统计显著性的低估（低于预期的假阳性率）。修正MVGC的两个版本都倾向于略微低估p值，从而导致略高于预期的假阳性率。在这里，没有显示标准MVGC模型的结果，该模型不包括特定的零项。因为数据是在零滞后时模拟的，MVGC没有建模，所以MVGC没有显示任何重要通道。这些结果支持了校正时间序列中自相关的必要性。

图5.包含时间自相关的模拟数据在所有方法上的(a)I型误差控制图和(b)ROC曲线。虚线表示方法的非稳健版本，而实线表示其稳健版本。

第二组模拟（图6）使用时间自相关和空间全局系统生理信号生成。正如预期的那样，AR相关性与Pearson相关系数相结合，导致ROC曲线中假阳性增加[图6(a)]和AUC降低[图6(b)]。当噪声独立于通道之间时，使用成对相关是有效的，但空间噪声引入了全局假阳性连接和模型显著性的大量过度报告。AR偏相关和修正的格兰杰因果关系仍然表现良好。在图5和图6中，模拟数据中没有添加运动伪影，因此稳健和非稳健方法在I型误差控制曲线和ROC曲线中表现相似。尽管使用稳健方法可能导致自由度损失，因为数据点的权重可能小于1，但时间点的数量仍然大得多（＞1000），因此自由度损失几乎可以忽略不计。

图6.包含时间自相关和全局信号的模拟数据在所有方法上的(a)I型误差控制图和(b)ROC曲线。虚线表示方法的非稳健版本，而实线表示其稳健版本。

在第三组模拟中，在包含时间自相关和全局信号的模拟数据中加入了以尖峰和位移伪影的头动。如图7所示，降低离群点权重的稳健方法比其非稳健方法的性能更好，具有稳健的修正格兰杰因果关系和AR偏相关，短距通道具有相似的AUC，但AR偏相关比修正MVGC具有更少的假阳性。值得注意的是，虽然稳健的修正格兰杰因果关系有效地减少了假阳性，但其仍然有远高于预期的假阳性率（0.05的预期率下，假阳性率为0.25）。

图7.包含时间自相关、全局信号和头动伪影信号的模拟数据在所有方法上的(a) I型误差控制图和(b)ROC曲线。

短距通道的有效性比较

本研究深入地探讨这些模型，并比较了仅用短距通道、长距通道或同时使用短距通道和长距通道对白化偏相关方法性能的影响。图8(a)显示，在I型误差控制图中，仅使用长距通道是最差的，可能会导致检验统计量（相关系数）显著性的高估，表明结果仍然保留了太多的全局连接。仅使用短距通道产生了最佳结果，尽管它略微低估了稳健估计器的显著性，并导致在估计中略微引入假阴性。最后，短距和长距通道的使用略差于单独使用短距通道，但优于仅使用长距通道。这可能是由于模型中的某些过拟合导致的。

图8(b)显示了附加运动伪影的模拟数据的比较。正如预期的那样，稳健的统计估计器大大优于非稳健的统计估计，仅用短距通道模型对这些伪影最敏感。与静止模拟一致，仅使用稳健的短距通道和短距和长距的稳健版本给出了最佳结果。在稳健模型中，仅使用长距通道仍然对测试统计的显著性有大量过度报告。

最后，图9显示了，使用修正的MVGC模型代替AR偏相关性进行了相同的分析。图9(a)显示了无运动模拟数据的结果。结果发现MVGC对模型中使用的通道不太敏感，短距、长距和甚至是两种通道模型都具有相似的性能。据观察，即使在静止数据模拟的情况下，MVGC的稳健版本对测试统计数据的显著性也有大量的过度报告。研究者认为，这是模型中关于估计有效自由度的未修正误差。特别是，由于测量变量y、x和/或Z中的时间点都可能是运动引起的异常值，因此难以在稳健估计器中定义统计异常值。该模型中使用的Tukey双平方权重假设离群值是独立的，但当前的模型没有有效地考虑多个通道的共享离群值时间点，因此自由度稍微高估。在图9(b)中，针对添加运动伪影的数据模拟比较MVGC模型，结果发现稳健方法优于非稳健方法。

建模滞后连接

修正后的MVGC同时考虑了零时滞和非零时滞。因此，当仅使用零时滞连接模拟数据时，AR偏相关性被确定为最佳。标准MVGC（非零时滞）在前几节中没有出现，是因为它在零非零时滞模拟中没有发现任何结果[图10(a)]。同样，当数据仅以非零时滞项的关系进行模拟时，AR偏相关模型失败，如图10(b)所示。修正后的MVGC是唯一一个在两种情况下都有效的模型。

在图11中，在添加了运动伪影的情况下，显示出了相同的模式，其中ROC曲线对于零时滞模拟的情况下的稳健AR偏相关和稳健修正MVGC模型是最佳的，并且对于非零时滞模拟下的稳健MVGC和稳健修正的MVGC是最佳的。在这两种情况下，修正后的MVGC都不如AR偏相关或MVGC好，因此如果知道存在哪种连接，则修正后的MVGC就不是首选。反之，这是两种情况下唯一可行的模型。

静息态fNIRS数据的结果

为了检验本研究的功能连接方法对实验静息态数据的适用性，研究者比较了Pearson相关、AR相关和AR偏相关的稳健版本，以在24名参与者样本上模拟通道之间的功能连接。图12显示了三种方法连接路径的组水平结果（q＜0.001，FDR校正）。结果表明，与AR相关和Pearson相关相比，使用AR偏相关导致更稀疏的网络，因为它校正了全局生理信号和自相关。使用Pearson相关性揭示了一个非常密集的连接网络，但大多数连接可能是假阳性的。虽然没有已知的地面真值或黄金标准的RS-fMRI网络可用于定量比较这些方法，但AR偏相关可能是三种方法中假阳性率最低的一种方法。

结论

复杂的数值模拟模型在理论上可以获取更多关于fNIRS时间通道之间关系的信息，但在实践中很难求解。尽管修正的格兰杰因果关系等方法可以对零时滞和时滞关系进行建模，但当连通性信息存在零时滞时，预白化偏相关性表现最佳，而当时间序列之间存在时滞关系时，格兰杰因果关系模型表现最佳。修正的MVGC在两种情况下都表现良好，但在低信噪比情况下，自由度较低，灵敏度可能更低。同样，带通滤波虽然可能会提高信噪比和灵敏度，但也可能导致自由度损失和自相关增加，从而抵消增益。因此，建议仔细考虑预处理流程中的步骤，以最大限度地提高灵敏度，同时将假阳性率降低到期望水平。

原文：Correction of global physiology in resting-state functional near-infrared spectroscopy.

DOI: 10.1117/1.NPh.9.3.035003

小伙伴们点个“在看”，加

(星标)关注茗创科技，将第一时间收到精彩内容推送哦～