導(dǎo)讀
意義：功能性近紅外光譜（fNIRS）數(shù)據(jù)的靜息態(tài)功能連接（RSFC）分析揭示了整個大腦的皮層連接和網(wǎng)絡(luò)汉形。fNIRS信號中的運動偽影和系統(tǒng)生理學(xué)提出了獨特的分析挑戰(zhàn)，目前已有研究探索了控制系統(tǒng)生理噪聲的方法弧哎。但是這些相同的方法在應(yīng)用于靜息態(tài)fNIRS（RS-fNIRS）數(shù)據(jù)時是否需要修改仍不清楚渠退。
目的：本研究系統(tǒng)地考察了幾種RSFC分析管道的敏感性和特異性庇忌，以確定校正RS-fNIRS數(shù)據(jù)中全局系統(tǒng)生理信號的最佳方法睛挚。
方法：使用數(shù)值模擬的RS-fNIRS數(shù)據(jù)，比較了幾種連通性分析管道的真陽性率和假陽性率溃蔫。使用受試者工作特征分析對其性能進行評分健提。管道包括偏相關(guān)和多變量格蘭杰因果關(guān)系，有/沒有短距測量伟叛，以及包括非傳統(tǒng)零滯后交叉項的修正多變量因果關(guān)系模型矩桂。研究者還考察了預(yù)白化和穩(wěn)健統(tǒng)計估計對性能的影響。
結(jié)果：與以往關(guān)于雙變量相關(guān)模型的研究一致痪伦，本研究結(jié)果表明穩(wěn)健統(tǒng)計和預(yù)白化是校正fNIRS時間序列中運動偽影和自相關(guān)的有效方法侄榴。此外，當(dāng)兩個感興趣的信號同步波動時网沾，使用從短距fNIRS通道提取的主成分作為偏相關(guān)模型的一部分進行預(yù)濾波癞蚕，在減少由于共有的系統(tǒng)生理學(xué)引起的偽相關(guān)方面最有效。然而辉哥，當(dāng)信號之間存在時間滯后時桦山，結(jié)合短距通道的多變量格蘭杰因果檢驗效果更好。由于不知道實驗數(shù)據(jù)中是否存在這樣的滯后醋旦，研究者提出了一種修正版本的格蘭杰因果關(guān)系恒水，其中包括非傳統(tǒng)的零滯后項作為折中解決方案。
結(jié)論：建議結(jié)合預(yù)白化饲齐、穩(wěn)健統(tǒng)計方法和處理管道中的偏相關(guān)性钉凌，以減少自相關(guān)、運動偽影和全局生理噪聲捂人，從而獲得具有統(tǒng)計有效性的RS-fNIRS連通性指標(biāo)御雕。未來的研究應(yīng)使用人類數(shù)據(jù)來驗證這些方法的有效性。
前言
自發(fā)或靜息態(tài)功能連接（RSFC）已成為表征人腦網(wǎng)絡(luò)結(jié)構(gòu)的重要工具滥搭。Biswal等人首先指出酸纲，功能性磁共振成像（fMRI）血氧水平依賴性（BOLD）信號中的同步自發(fā)波動被組織成不同的網(wǎng)絡(luò)，這些網(wǎng)絡(luò)在被試之間保持不變瑟匆。25年后闽坡，這些結(jié)果在數(shù)百項研究中被復(fù)制，并擴展到人類大腦功能、發(fā)育和疾病的各個方面疾嗅。靜息態(tài)網(wǎng)絡(luò)也使用其他神經(jīng)成像方法進行了觀察外厂，包括EEG和MEG，以及fNIRS宪迟。靜息態(tài)范式的一個主要優(yōu)點是不需要執(zhí)行任務(wù)，這在嬰兒交惯、老年人次泽、患者等人群的研究中特別有用。
目前已經(jīng)開發(fā)了幾種分析方法來量化靜息態(tài)連通性席爽。對于fMRI和fNIRS意荤，這些方法通常關(guān)注大腦區(qū)域之間的零滯后關(guān)系。零滯后關(guān)系是指兩個信號同時發(fā)生信號波動（在測量采樣時間內(nèi)）的關(guān)系只锻。估計功能連通性的最常用方法可能是計算兩個腦區(qū)時間序列之間的Pearson相關(guān)系數(shù)玖像。有幾種方法可以在時域中表征RSFC，包括（i）基于種子的連接齐饮，（ii）成對連接捐寥，以及（iii）獨立成分分析（ICA）。在基于種子的連接中祖驱，選擇一個種子區(qū)域握恳，并評估其與所有其他腦區(qū)域的連通性。在成對連接中捺僻，計算所有可能通道對之間的相關(guān)性乡洼。ICA將fNIRS通道數(shù)據(jù)的協(xié)方差分解為幾個統(tǒng)計上獨立的時空成分，其中一些成分表示每個被試的大腦功能網(wǎng)絡(luò)匕坯，或者表示跨被試連接的功能網(wǎng)絡(luò)束昵。
格蘭杰因果關(guān)系是“有效連接”方法的一個例子，該方法依賴于fNIRS通道在時域上的滯后關(guān)系葛峻。另一種常用于模擬fNIRS中功能連通性的方法是小波變換相干（WTC）锹雏，該方法在時頻域中模擬通道之間的關(guān)系。WTC依賴于連續(xù)小波變換（CWT）术奖，CWT將時間序列數(shù)據(jù)轉(zhuǎn)換為時頻域逼侦。WTC估計為時間序列在時頻域中的互相關(guān)。
功能性近紅外光譜
fNIRS是一種光學(xué)成像方法腰耙，使用近紅外波長范圍（650至950nm）的光來估計氧合血紅蛋白（HbO）和脫氧血紅蛋白（HbR）濃度的時間變化榛丢。fNIRS測量大腦的血流動力學(xué)變化，反映血流量挺庞、體積和血氧合的變化晰赞，因此，它是一種間接測量神經(jīng)活動的方法，就像fMRI BOLD信號一樣掖鱼。與fMRI相比然走，fNIRS具有一些優(yōu)勢，包括成本更低戏挡，對頭部運動的彈性更大芍瑞，以及更高的數(shù)據(jù)采樣率。fNIRS也可以便攜褐墅，能夠執(zhí)行更靈活或生態(tài)有效性的任務(wù)拆檬。然而，fNIRS的缺點是其空間分辨率較低妥凳，深入大腦的深度有限竟贯，相對于大腦本身，對頭部表層的敏感性增加逝钥。對于fNIRS測量來說屑那，后一個問題尤其令人困惑，因為頭皮的系統(tǒng)生理信號艘款，包括心臟持际、呼吸和血壓波動，通常會導(dǎo)致空間上的全局噪聲哗咆，其振幅通常大于感興趣的潛在大腦信號选酗。因此，在使用fNIRS的RSFC實驗中岳枷，將表層生理噪聲和真正的潛在神經(jīng)相關(guān)腦網(wǎng)絡(luò)進行分離是具有挑戰(zhàn)性的芒填。
Obrig等人首次證明了fNIRS在捕獲極低頻和低頻振蕩方面的可行性，這些振蕩有助于RSFC的研究空繁。之后殿衰，幾項用于考察RSFC可行性的早期研究使用fNIRS和擴散光學(xué)層析成像（DOT）發(fā)現(xiàn)了可靠的和預(yù)期的視覺、感覺運動和語言網(wǎng)絡(luò)盛泡。從那時起闷祥，研究人員使用靜息態(tài)fNIRS（RS-fNIRS）和DOT研究嬰兒發(fā)育、前額葉皮層的性別差異傲诵、早產(chǎn)凯砍、語言網(wǎng)絡(luò)、自閉癥譜系障礙拴竹、情感障礙和衰老悟衩。檢查RSFC的重測信度發(fā)現(xiàn)，在從一小時到一周的多個sessions中栓拜，個體和群體水平的RSFC圖之間存在良好到極好的組內(nèi)相關(guān)性座泳。
RS-fNIRS信號的統(tǒng)計特性
鑒于RS-fNIRS的信號和噪聲特性惠昔，與RS-fMRI相比，RS-fNIRS具有一些獨特的挑戰(zhàn)挑势。由于fNIRS信號的采集方式镇防，它對頭皮的非神經(jīng)腦外信號變化很敏感。RS-fNIRS信號由三部分組成：（i）低頻神經(jīng)振蕩（感興趣信號）潮饱，（ii）源自大腦的腦內(nèi)生理噪聲来氧，以及（iii）源自頭皮和非腦組織的腦外生理噪聲。鑒于靜息態(tài)范式的性質(zhì)香拉，所有信號都是非誘發(fā)的啦扬。因此，fNIRS信號同時受到腦外信號變化和非神經(jīng)腦信號變化的影響缕溉，包括心臟（約1-1.2Hz）考传、呼吸（0.3-0.6Hz）和血壓/Mayer波（0.1Hz）波動引起的系統(tǒng)生理噪聲吃型。生理噪聲對RS-fNIRS信號有兩種影響：①誘導(dǎo)時間自相關(guān)和②增加大腦通道之間的空間協(xié)方差证鸥。
先前的研究表明，高采樣率和系統(tǒng)生理信號(心率勤晚、呼吸和血壓)枉层，結(jié)合潛在的血流動力學(xué)反應(yīng)，可以共同導(dǎo)致時序相關(guān)赐写，從而在RS-fNIRS時間序列中會產(chǎn)生自相關(guān)鸟蜡。去除時間序列中的自相關(guān)以獲得統(tǒng)計有效的關(guān)系估計，可以使用預(yù)白化的方法挺邀。預(yù)白化是通過去除時間序列中的序列相關(guān)來實現(xiàn)的揉忘。自回歸(AR)模型可用于數(shù)據(jù)預(yù)白化和去除序列相關(guān)。具體來說端铛，對于fNIRS泣矛，使用預(yù)白化需要使用p階AR模型來模擬當(dāng)前時間點作為前一個時間點的函數(shù)。AR模型中去除時間相關(guān)性后的殘差項被用來計算通道之間的Pearson相關(guān)系數(shù)禾蚕，而不是使用自相關(guān)通道數(shù)據(jù)您朽。當(dāng)將預(yù)白化納入到預(yù)處理管道中用于功能連接估計時，可以減少I型錯誤换淆。盡管對fNIRS的時間自相關(guān)進行了校正哗总，但如果沒有考慮由于系統(tǒng)生理機能導(dǎo)致的空間協(xié)方差增加，這也可能會導(dǎo)致FDR的增加倍试。
由于系統(tǒng)生理信號被認為在整個大腦中是相對均勻的讯屈，所以由所有通道共有的系統(tǒng)波動引起的信號變化可以增加通道中的空間協(xié)方差。圖1顯示县习，不校正fNIRS數(shù)據(jù)中增加的空間協(xié)方差可能會導(dǎo)致I型誤差增加和特異性降低耻煤，因為基于地面真值的連接路徑遠比預(yù)期的多具壮。
方法
為了表征每種分析方法的敏感性和特異性怜跑，本研究根據(jù)具有已知地面真值的數(shù)值模擬數(shù)據(jù)，探索了模型在幾種模擬情況中的性能吠勘，從簡單正態(tài)分布隨機噪聲（其中本研究模型的所有統(tǒng)計假設(shè)都有效）到具有時間自相關(guān)噪聲的生理信號性芬，再到全局空間生理噪聲。最后剧防，加入運動偽影植锉。每個數(shù)據(jù)模擬都代表了一個更具挑戰(zhàn)性的統(tǒng)計問題，通過模擬數(shù)據(jù)中的真實和零相關(guān)性峭拘，使用受試者工作特征（ROC）分析來檢驗這個問題俊庇。所有方法均在MATLAB（R2019b）的NIRS AnalyzIR工具箱中實現(xiàn)。
數(shù)據(jù)模擬
采用具有指定空間協(xié)方差矩陣的多變量正態(tài)分布方法鸡挠，對由16個長距通道（27mm）和9個短距通道（7.1mm）組成的探測器進行了fNIRS數(shù)據(jù)模擬辉饱。模擬了近紅外波長為690和830nm的數(shù)據(jù)。這種已知的空間協(xié)方差在數(shù)據(jù)中產(chǎn)生潛在的“真實”相關(guān)性拣展，在此基礎(chǔ)上評估各種算法彭沼。這些真實相關(guān)性是通過通道之間的零滯后（無時間偏移）或時間滯后（數(shù)學(xué)因果模型）產(chǎn)生的。僅針對長距通道在通道空間中模擬目標(biāo)真實相關(guān)性备埃。為了給模型添加額外的噪聲姓惑，使用半無限均勻平板模型來計算光學(xué)正向模型。使用光學(xué)正向模型將這些體素投影到通道空間按脚，并添加到基本真值相關(guān)信號于毙。
對于每個模擬，以4Hz采樣率生成300 s的靜息態(tài)數(shù)據(jù)乘寒。為了生成ROC曲線望众，每個測試總共生成了200個數(shù)據(jù)集。在ROC分析中伞辛，在探測器中任何兩個通道之間的所有可能連接中烂翰，模擬了10%左右的可能連接路徑的真實連接，確保全局系統(tǒng)噪聲是主要的空間噪聲特征蚤氏。因此甘耿，用于生成模擬數(shù)據(jù)的協(xié)方差矩陣中只有10%的非對角元素為非零。由于時間自相關(guān)竿滨、全局系統(tǒng)生理信號或運動偽影佳恬，其余90%路徑之間的任何顯著連接都是假陽性捏境。為了生成真和假陽性估計，從總鄰接矩陣中隨機抽取相等數(shù)量的真和空連接毁葱。
fNIRS實驗數(shù)據(jù)
采集
研究者還測試了本研究的功能連接方法對實驗靜息態(tài)數(shù)據(jù)的有效性和適用性垫言。使用商業(yè)NIRScout-2（NIRx Medizintechnik GmbH，柏林倾剿，德國）連續(xù)波fNIRS系統(tǒng)記錄NIRS數(shù)據(jù)筷频。共有50個通道（42個長距通道，8個短距通道）分布在雙側(cè)額葉和感覺運動皮層（見圖2）前痘。如圖2所示凛捏，長距通道由放置在頭皮上的16個源光電二極管(橙色圓圈)和13個探測器(藍色矩形)組成。8個探測器（綠色菱形）放置在頭皮的8個短距通道上芹缔。藍色實線表示長距通道坯癣，綠色虛線表示短距通道。以7.8125Hz的采樣率記錄兩個波長（760nm和850nm）的數(shù)據(jù)最欠。定位頭帽后示罗，使用NIRx Aurora軟件優(yōu)化信號質(zhì)量。
參與者和任務(wù)
約24名健康被試參與了實驗（9名男性和15名女性鹉勒；平均年齡27.4歲帆锋，SD=8.7歲吵取；23名右利手）。告知被試實驗情況锯厢，并獲得書面同意皮官。這項研究得到了匹茲堡大學(xué)機構(gòu)審查委員會的批準(zhǔn)。每個被試進行一次靜息態(tài)掃描和其他任務(wù)掃描。指示被試在靜息態(tài)掃描期間盡量減少身體運動，保持坐姿放松5分鐘鞠绰，而不使用任何腦力勞動绍些。
數(shù)據(jù)預(yù)處理
將模擬的fNIRS信號強度轉(zhuǎn)換為光密度的變化，然后使用修正的Beer-Lambert定律（MBLL）將光密度單位轉(zhuǎn)換為HbO和HbR濃度瑰煎，兩種波長的微分路徑長度因子為6，部分體積校正均為60。除了對實驗fNIRS時間序列進行降采樣外馍佑，以與模擬數(shù)據(jù)類似的方式對數(shù)據(jù)進行預(yù)處理。使用MBLL將原始信號強度轉(zhuǎn)換為HbO和HbR濃度后梨水，將實驗靜息態(tài)fNIRS數(shù)據(jù)從7.8125Hz降采樣至4Hz拭荤，以匹配模擬數(shù)據(jù)的采樣率。在本文的其余部分疫诽，分析的重點在于HbO濃度的變化舅世，結(jié)果也應(yīng)適用于HbR旦委。
圖3顯示了三種情況下模擬數(shù)據(jù)的fNIRS時間序列示例，（i）僅時間自相關(guān)雏亚，（ii）具有時間自相關(guān)和全局系統(tǒng)生理信號缨硝，以及（iii）具有時間自相關(guān)、全局系統(tǒng)生理和運動偽影罢低。圖3(左)顯示了在使用MBLL將原始模擬fNIRS信號轉(zhuǎn)換為HbO信號變化后追葡，三種情況下的兩個通道的時間軌跡。圖3(右)顯示了通道的相關(guān)結(jié)構(gòu)奕短。正如預(yù)期的那樣宜肉，隨著全局系統(tǒng)生理信號的引入，這些通道中的協(xié)方差增加翎碑。隨著在模擬數(shù)據(jù)中引入運動偽影谬返，隨著通道之間的相關(guān)值接近1，協(xié)方差進一步急劇增加日杈。
預(yù)白化
預(yù)白化通常使用AR模型來去除時間序列中的自相關(guān)，并對信號的頻率內(nèi)容進行白化莉擒。預(yù)白化大大降低了連接模型中的FDR酿炸。這里，使用AR模型實現(xiàn)預(yù)白化涨冀，如下所示：
fNIRS時間進程可表示為：
其中εt∈N(0填硕，σ2)，yt是fNIRS的時間序列鹿鳖，ai是AR參數(shù)扁眯，p是模型階數(shù)，由貝葉斯信息準(zhǔn)則（BIC）確定翅帜，εt是新息項姻檀，用于計算連通性指標(biāo)，而不是時間序列本身涝滴。圖4顯示了三組模擬數(shù)據(jù)預(yù)白化前后代表性通道的（a）功率譜密度和（b）自相關(guān)函數(shù)圖绣版，包括：（i）時間自相關(guān)，（ii）時間自相關(guān)和全局系統(tǒng)生理信號歼疮，（iii）時間自相關(guān)杂抽、全局系統(tǒng)生理和運動偽影。正如預(yù)期的那樣腋妙，預(yù)白化步驟使功率譜密度圖變得平坦默怨，并降低了fNIRS時間序列中的時間自相關(guān)。
穩(wěn)健方法
在對時間序列進行預(yù)白化后匙睹，新息項包含在該時間點添加的信息愚屁。由于運動偽影被視為異常值，它們可能會影響正常性假設(shè)痕檬。通常霎槐，這些高影響點會對連通性指標(biāo)的估計產(chǎn)生重大影響。如前所述梦谜，運動偽影對連接性指標(biāo)的影響不同丘跌，這取決于運動是跨多個通道共享還是僅限于幾個通道⊙渥跨多個通道共享的運動損壞時間點可以通過增加它們之間的協(xié)方差來影響功能連接闭树。因此，為了識別統(tǒng)計異常值荒澡，研究者使用了之前在Santosa等人研究中引入的聯(lián)合加權(quán)函數(shù)的多變量版本报辱。聯(lián)合加權(quán)函數(shù)可以通過計算時間序列的幾何長度來識別異常值，并根據(jù)它們偏離平均值的方式對它們進行適當(dāng)?shù)慕禉?quán)单山，將運動損壞的異常值分配較低的權(quán)重碍现，使它們對相關(guān)和回歸方法的影響較低。
根據(jù)所有通道（包括短距通道）的新息項計算幾何長度函數(shù)rt米奸，如下所示：
其中εi昼接，t表示時間點t的第i個通道的新息項，c是通道總數(shù)悴晰。
根據(jù)離群值偏離平均值的程度對離群值進行加權(quán)的加權(quán)函數(shù)(S)慢睡，可使用Tukey雙平方函數(shù)的平方根，對所有時間點進行如下計算：
這里膨疏，σ是分布的標(biāo)準(zhǔn)偏差一睁，根據(jù)中值絕對偏差（MAD）和正態(tài)分布標(biāo)準(zhǔn)偏差（σ）之間的關(guān)系估計(σ≈ 1.4826*MAD）钻弄。使用MAD是因為它是更穩(wěn)健的分散統(tǒng)計佃却。常數(shù)κ的值為4.685，這對于沒有任何異常值的標(biāo)準(zhǔn)正態(tài)分布提供了95%的效率窘俺。具有對角元素的加權(quán)矩陣S表示所有時間點的加權(quán)函數(shù)可與所有時間序列的新息項相乘饲帅，從而降低正態(tài)分布統(tǒng)計異常點的權(quán)值。
其中瘤泪，Yc是所有通道重加權(quán)后的新息項矩陣灶泵，Ec是預(yù)白化后獲得的新息項矩陣。
因此对途，為了減少運動偽影的影響赦邻，對于回歸模型，不使用線性模型和普通最小二乘估計參數(shù)实檀，而是使用穩(wěn)健統(tǒng)計惶洲，使用加權(quán)最小二乘（WLS）來降低從學(xué)生化殘差計算出的異常值的權(quán)重按声。權(quán)重是使用Tukey的雙平方函數(shù)計算的。然后將計算出的權(quán)重應(yīng)用于原始數(shù)據(jù)恬吕，并重新估計殘差签则。重復(fù)該過程直到收斂。穩(wěn)健回歸在內(nèi)置的MATLAB函數(shù)“robustfit”中實現(xiàn)铐料。MATLAB函數(shù)robustfit使用50次迭代或回歸系數(shù)值無可檢測數(shù)值變化的收斂標(biāo)準(zhǔn)來評估系數(shù)的收斂性渐裂。由于多變量格蘭杰因果關(guān)系（MVGC）和修正的MVGC在50次迭代中不收斂，研究者將最大迭代限制更改為150次迭代钠惩，以確逼饬梗回歸參數(shù)收斂。
使用WLS回歸估計的參數(shù)可用如下公式獲得：
這里篓跛，是參數(shù)估計扛拨，S是加權(quán)矩陣。
用于減少全局信號變化的基于PCA的回歸
如果探針廣泛分布在大腦中举塔，那么跨多個通道共享的大信號方差可能主要反映生理噪聲绑警。因此，PCA的前幾個成分可用于識別由于生理噪聲引起的全局信號變化央渣，并從通道數(shù)據(jù)中回歸其影響计盒。數(shù)據(jù)協(xié)方差矩陣可以使用奇異值分解法進行分解
其中C是給定C=YT·Y的數(shù)據(jù)協(xié)方差矩陣，U是包含左奇異向量的矩陣芽丹”逼簦∑是包含奇異值的對角矩陣，V是包含右奇異向量的矩陣拔第。通常咕村，使用解釋最大方差的前幾個成分。研究者使用最小數(shù)量的主成分來解釋至少90%的通道協(xié)方差蚊俺。對長距通道懈涛、短距通道或兩者進行協(xié)方差分解可以提取全局信號。這三種選項都在測試中進行了檢查泳猬。經(jīng)PCA回歸后的通道數(shù)據(jù)如下：
這里批钠，y是新息項或通道的原始數(shù)據(jù)，取決于分析通道得封。Z是包含前幾個成分作為列的矩陣埋心。
功能連接指標(biāo)
本研究考察了幾種估計功能連接的方法。在這里忙上，重點關(guān)注時域中的方法拷呆，包括相關(guān)性、偏相關(guān)和格蘭杰因果關(guān)系的幾種變式。時頻方法茬斧，如小波相干箫柳，由于兩個原因沒有在本研究中進行檢查。首先啥供，在研究者之前的工作中發(fā)現(xiàn)悯恍，在雙變量連接性分析中，小波相干優(yōu)于ROC分析中的標(biāo)準(zhǔn)Pearson相關(guān)伙狐，但性能不如預(yù)白化的穩(wěn)健相關(guān)模型涮毫。除非標(biāo)準(zhǔn)小波相干模型被修改為包括高階AR項和穩(wěn)健（異常值拒絕）方法贷屎，否則小波相干也具有較高的假陽性率和不受控的I型誤差罢防。經(jīng)過這些修改發(fā)現(xiàn)，與預(yù)白化穩(wěn)健模型相比唉侄，小波相干沒有任何優(yōu)勢咒吐。第二，小波相干的多變量或偏相關(guān)版本尚未描述属划，需要開發(fā)一個新的統(tǒng)計模型恬叹。
基于相關(guān)性的分析模型
Pearson相關(guān)系數(shù)。兩個時間序列x和y之間的Pearson相關(guān)系數(shù)計算如下：
其中σx同眯，y是時間序列x和y的協(xié)方差绽昼，σx和σy分別是x和y標(biāo)準(zhǔn)差。標(biāo)準(zhǔn)Pearson相關(guān)性并不是一個穩(wěn)健的統(tǒng)計估計须蜗，Santosa等人在fNIRS的背景下提出并描述了若干穩(wěn)健變量硅确。穩(wěn)健相關(guān)性計算為回歸模型y對x和x對y的標(biāo)準(zhǔn)回歸系數(shù)的幾何平均值，以減少x或y中的異常值明肮，但這仍然會影響相關(guān)估計菱农。步驟如下：
基于上述回歸模型中的學(xué)生化殘差估計加權(quán)矩陣，加權(quán)矩陣Sx→y在上述回歸模型的兩側(cè)相乘柿估，并且更新參數(shù)估計循未。
同樣，可以在y上擬合x官份，
加權(quán)矩陣Sy→x根據(jù)學(xué)生化殘差估計只厘，以降低x中的異常值。Sy→x在上述回歸模型的兩側(cè)相乘舅巷，并更新參數(shù)估計。
兩個回歸模型交替求解河咽，直到參數(shù)βy→x和βx→y收斂钠右。使用以下公式估計最終穩(wěn)健相關(guān)系數(shù)：
穩(wěn)健相關(guān)系數(shù)估計的符號由回歸系數(shù)的符號確定。有效自由度為
忘蟹。通過降低數(shù)據(jù)點的權(quán)重飒房，特別是對異常點的權(quán)重搁凸，降低了有效自由度。
對于標(biāo)準(zhǔn)和穩(wěn)健相關(guān)估計狠毯，在給定零假設(shè)（p值）的情況下护糖，觀察數(shù)據(jù)的概率可以通過使用T統(tǒng)計量估計
其中df為標(biāo)準(zhǔn)Pearson相關(guān)的有效自由度，值為n?2嚼松，r為Pearson相關(guān)系數(shù)嫡良。這里n是時間點的數(shù)量。
自回歸相關(guān)献酗。AR相關(guān)性是基于預(yù)白化時間序列計算的Pearson相關(guān)系數(shù)寝受。如Santosa等人所述，為每個數(shù)據(jù)通道估計p階AR模型罕偎，并用于預(yù)白信號很澄。然后在相關(guān)分析中使用AR模型中的新息項（例如，在每個時間步長添加的新信息）颜及。研究者使用的最大模型階數(shù)為40（10×采樣率）甩苛，高于模擬fNIRS數(shù)據(jù)時使用的AR濾波器。根據(jù)經(jīng)驗觀察俏站，在實驗靜息態(tài)fNIRS數(shù)據(jù)中浪藻，模型階數(shù)通常＜10·Fs（此處Fs為采樣頻率）。
AR偏相關(guān)乾翔。偏相關(guān)涉及一個兩步過程爱葵，其中首先將其他通道從感興趣的數(shù)據(jù)中投影出來，然后對殘差執(zhí)行相關(guān)反浓。使用基于PCA回歸的投影來部分消除其他通道的影響萌丈。PCA回歸用于避免模型的共線不穩(wěn)定性。使用最小數(shù)量的特征向量來解釋至少90%的空間協(xié)方差雷则。對于兩個通道i和j辆雾，如果預(yù)白化后的新息項由εi和εj表示，則白化PCA濾波時間序列計算為
現(xiàn)在月劈，當(dāng)從通道中去除全局信號時度迂，AR偏相關(guān)可以計算為兩個殘差之間的Pearson相關(guān)性〔戮荆基于保留的主成分數(shù)量進一步減少自由度惭墓。該模型的輸入數(shù)據(jù)是預(yù)白化新息項時間序列，而不是通道數(shù)據(jù)而姐，因此該方法被稱為AR偏相關(guān)腊凶。
基于格蘭杰因果關(guān)系的分析模型
多變量格蘭杰因果關(guān)系。在控制其他時間序列的影響后，MVGC檢驗兩個時間序列之間的滯后關(guān)系钧萍。這本質(zhì)上是兩個回歸模型的模型擬合檢驗褐缠。可以使用多變量AR模型（MVAR）控制通道Z后檢驗x“格蘭杰原因”是否為y风瘦。對于通道時間序列队魏，可以將當(dāng)前時間點yt建模為限制性模型中y和其他通道Z之前滯后的函數(shù)：
（限制性模型）
在非限制性模型中，添加x的滯后項万搔，以模擬y中的當(dāng)前時間點：
（非限制性模型）
這里胡桨，p是模型中包含的滯后項的數(shù)量，使用BIC的擬合測量確定蟹略。為了檢驗因果關(guān)系登失，進行了以下統(tǒng)計檢驗，與沒有x（限制性模型）的模型相比挖炬，在x（非限制性模型）中添加額外滯后項是否提高了yt的可預(yù)測性揽浙。如果x的過去包含關(guān)于當(dāng)前時間點yt的預(yù)測信息，則即使在控制丟失的自由度之后意敛，yt預(yù)測中的誤差也會得到改善馅巷。這在統(tǒng)計上可以評估為：
零假設(shè)：x不是y的格蘭杰原因
H0: β1=β2=β3=……βp=0
備擇假設(shè)：x是y的格蘭杰原因。
HA：至少有一個βi為非零
因為有嵌套模型草姻，所以可以使用嵌套F統(tǒng)計量來檢驗格蘭杰因果關(guān)系
這里钓猬，SSER和SSEU分別是限制性和非限制性模型的平方誤差之和，n是時間點的數(shù)量撩独。因此敞曹，可以在自由度p和[n?(2p+1)]處計算F統(tǒng)計量的p值。格蘭杰因果指標(biāo)G也可以用于評估滯后關(guān)系的強度综膀，如下所示：
在雙變量格蘭杰因果關(guān)系上使用MVGC的目的是控制包括系統(tǒng)生理信號在內(nèi)的第三變量的影響澳迫。在MVGC模型中添加多個通道的滯后項，容易使回歸模型無法求解剧劝，因為變量的數(shù)量將超過時間點的數(shù)量橄登。因此，與使用通道數(shù)據(jù)本身不同讥此，首先執(zhí)行PCA分析拢锹，并像之前使用偏相關(guān)方法一樣，包含解釋至少90%方差的最少成分萄喳。然后卒稳，在限制性和非限制性模型中包括成分和時滯成分，而不是通道數(shù)據(jù)取胎。
修正的多變量格蘭杰因果關(guān)系展哭。傳統(tǒng)的MVGC僅包括方程右側(cè)的y湃窍、x和Z的時間滯后歷史闻蛀。因此匪傍，MVGC僅對非零滯后關(guān)系建模。重要的是要認識到零滯后與采樣率有關(guān)觉痛，在模擬fNIRS數(shù)據(jù)的4Hz采樣時役衡，一個采樣（250ms）內(nèi)發(fā)生的任何同步信號變化在數(shù)學(xué)上都是零滯后。雖然大腦中的神經(jīng)信號傳遞相當(dāng)快（100ms）薪棒，但我們不希望假設(shè)所有的相關(guān)性都是零滯后手蝎，因而本研究引入了一個改進的MVGC模型，以建模通道之間的滯后和零滯后關(guān)系俐芯。在修正的MVGC中棵介，將Z的零滯后項包括在限制性和非限制性模型中，以及通道x的零滯后項包括在非限制性模型中吧史。因此邮辽，限制性模型可以寫成：
（限制性模型）
類似地，非限制性模型可以寫成：
（非限制性模型）
上述方程還包括Z和x的零滯后項∶秤現(xiàn)在吨述，為了檢驗格蘭杰因果關(guān)系的統(tǒng)計顯著性，檢驗了如下零假設(shè)和備擇假設(shè)：
H0: β1=β2=β3=……βp=0
HA：至少有一個βi為非零
統(tǒng)計顯著性采用嵌套F統(tǒng)計進行評估钞脂，格蘭杰因果關(guān)系的計算與傳統(tǒng)MVGC類似揣云。該修正模型的一個警告是，假設(shè)不檢驗y和x之間的連接是否是由于零或非零滯后項冰啃。我們可以通過修改非限制性模型以包含x的非零滯后項來進一步檢驗這一點邓夕，使得非限制性模型和限制性模型之間的唯一區(qū)別是在x中包含零滯后項。
穩(wěn)健的多變量格蘭杰因果關(guān)系阎毅。最后焚刚，研究者開發(fā)了MVGC和修正的MVGC模型的穩(wěn)健統(tǒng)計版本。該多變量模型中離群值的定義具有挑戰(zhàn)性净薛，因為如雙變量相關(guān)模型汪榔，離群值可能存在于y、x和/或Z的時間過程中肃拜。在雙變量相關(guān)中痴腌，使用離群值權(quán)重的并集（例如，y或x中的離群值）燃领。然而結(jié)果發(fā)現(xiàn)士聪，擴展到多個通道的數(shù)據(jù)通常會導(dǎo)致太多的數(shù)據(jù)被向下加權(quán)，因為任何一個通道中的異常點會導(dǎo)致所有通道的時間點被向下加權(quán)猛蔽。特別是剥悟，MVGC中的時滯要求包含該列矩陣的所有滯后列也被向下加權(quán)灵寺。
為了降低上述回歸模型中使用多個項的穩(wěn)健回歸方法的復(fù)雜性，研究者使用新息項区岗，而不是單獨使用通道數(shù)據(jù)略板。穩(wěn)健MVGC的過程如下所述。首先慈缔，對通道數(shù)據(jù)（短距或長距）進行預(yù)白化叮称，并計算新息項。使用從新息項計算出的加權(quán)矩陣對數(shù)據(jù)進行預(yù)加權(quán)藐鹤。經(jīng)過預(yù)白化和預(yù)加權(quán)后瓤檐，其他通道（Z）的新息項包括在限制性模型的回歸中，如下所示：
（限制性模型）
這里娱节，εty是在預(yù)白化和預(yù)加權(quán)后的時間點t處通道y的新息項挠蛉。εjz,t?i是時間點t-i處第j個成分的主成分權(quán)值。nPC表示使用的主成分的數(shù)量肄满。
同樣谴古，新息項也用于非限制性模型，而不是通道數(shù)據(jù)本身：
（非限制性模型）悄窃。
代碼實現(xiàn)
所有方法均在MATLAB中實現(xiàn)讥电。利用MATLAB內(nèi)置函數(shù)和NIRS Brain AnalyzIR工具箱對相關(guān)方法進行預(yù)處理和實現(xiàn)。使用自定義的MATLAB腳本實現(xiàn)了偏相關(guān)轧抗、格蘭杰因果關(guān)系和穩(wěn)健格蘭杰因果等方法恩敌。本文中使用的工具箱和代碼可在線訪問：https://github.com/pradlanka/rsfc-fnirs獲得。
結(jié)果
處理方法的比較
如預(yù)期那樣横媚，隨著模擬的進展和噪聲的增加纠炮，標(biāo)準(zhǔn)Pearson相關(guān)等基本模型很快就難以為繼了。例如灯蝴，Pearson相關(guān)假設(shè)獨立的測量恢口，這與加入的時間自相關(guān)噪聲相違背。如圖5(b)所示穷躁，與預(yù)白化和MVGC模型相比耕肩，不包括AR（預(yù)白化）濾波器的Pearson和穩(wěn)健Pearson相關(guān)模型在ROC曲線上的性能顯著更低。幾乎所有其他方法的AUC接近1问潭。此外猿诸，如圖5(a)所示，Pearson和穩(wěn)健Pearson模型具有較大的假陽性率狡忙。其他分析模型的I型誤差圖更接近理想（對角線）梳虽。AR相關(guān)和穩(wěn)健AR相關(guān)都具有最接近理想的I型誤差控制。AR偏相關(guān)模型略微高估了測試統(tǒng)計的顯著性灾茁，而穩(wěn)健的AR偏相關(guān)過校正導(dǎo)致了對測試統(tǒng)計顯著性的低估（低于預(yù)期的假陽性率）窜觉。修正MVGC的兩個版本都傾向于略微低估p值谷炸，從而導(dǎo)致略高于預(yù)期的假陽性率。在這里禀挫，沒有顯示標(biāo)準(zhǔn)MVGC模型的結(jié)果旬陡，該模型不包括特定的零項。因為數(shù)據(jù)是在零滯后時模擬的特咆，MVGC沒有建模季惩，所以MVGC沒有顯示任何重要通道录粱。這些結(jié)果支持了校正時間序列中自相關(guān)的必要性腻格。
圖5.包含時間自相關(guān)的模擬數(shù)據(jù)在所有方法上的(a)I型誤差控制圖和(b)ROC曲線。虛線表示方法的非穩(wěn)健版本啥繁，而實線表示其穩(wěn)健版本菜职。
第二組模擬（圖6）使用時間自相關(guān)和空間全局系統(tǒng)生理信號生成。正如預(yù)期的那樣旗闽，AR相關(guān)性與Pearson相關(guān)系數(shù)相結(jié)合酬核，導(dǎo)致ROC曲線中假陽性增加[圖6(a)]和AUC降低[圖6(b)]。當(dāng)噪聲獨立于通道之間時适室，使用成對相關(guān)是有效的嫡意，但空間噪聲引入了全局假陽性連接和模型顯著性的大量過度報告。AR偏相關(guān)和修正的格蘭杰因果關(guān)系仍然表現(xiàn)良好捣辆。在圖5和圖6中蔬螟，模擬數(shù)據(jù)中沒有添加運動偽影，因此穩(wěn)健和非穩(wěn)健方法在I型誤差控制曲線和ROC曲線中表現(xiàn)相似汽畴。盡管使用穩(wěn)健方法可能導(dǎo)致自由度損失旧巾，因為數(shù)據(jù)點的權(quán)重可能小于1，但時間點的數(shù)量仍然大得多（＞1000）忍些，因此自由度損失幾乎可以忽略不計鲁猩。
圖6.包含時間自相關(guān)和全局信號的模擬數(shù)據(jù)在所有方法上的(a)I型誤差控制圖和(b)ROC曲線。虛線表示方法的非穩(wěn)健版本罢坝，而實線表示其穩(wěn)健版本廓握。
在第三組模擬中，在包含時間自相關(guān)和全局信號的模擬數(shù)據(jù)中加入了以尖峰和位移偽影的頭動嘁酿。如圖7所示隙券，降低離群點權(quán)重的穩(wěn)健方法比其非穩(wěn)健方法的性能更好，具有穩(wěn)健的修正格蘭杰因果關(guān)系和AR偏相關(guān)痹仙，短距通道具有相似的AUC是尔，但AR偏相關(guān)比修正MVGC具有更少的假陽性。值得注意的是开仰，雖然穩(wěn)健的修正格蘭杰因果關(guān)系有效地減少了假陽性拟枚，但其仍然有遠高于預(yù)期的假陽性率（0.05的預(yù)期率下薪铜，假陽性率為0.25）。
圖7.包含時間自相關(guān)恩溅、全局信號和頭動偽影信號的模擬數(shù)據(jù)在所有方法上的(a) I型誤差控制圖和(b)ROC曲線隔箍。
短距通道的有效性比較
本研究深入地探討這些模型，并比較了僅用短距通道脚乡、長距通道或同時使用短距通道和長距通道對白化偏相關(guān)方法性能的影響蜒滩。圖8(a)顯示，在I型誤差控制圖中奶稠，僅使用長距通道是最差的俯艰，可能會導(dǎo)致檢驗統(tǒng)計量（相關(guān)系數(shù)）顯著性的高估，表明結(jié)果仍然保留了太多的全局連接锌订。僅使用短距通道產(chǎn)生了最佳結(jié)果竹握，盡管它略微低估了穩(wěn)健估計器的顯著性，并導(dǎo)致在估計中略微引入假陰性辆飘。最后啦辐，短距和長距通道的使用略差于單獨使用短距通道，但優(yōu)于僅使用長距通道蜈项。這可能是由于模型中的某些過擬合導(dǎo)致的芹关。
圖8(b)顯示了附加運動偽影的模擬數(shù)據(jù)的比較。正如預(yù)期的那樣紧卒，穩(wěn)健的統(tǒng)計估計器大大優(yōu)于非穩(wěn)健的統(tǒng)計估計侥衬，僅用短距通道模型對這些偽影最敏感。與靜止模擬一致常侦，僅使用穩(wěn)健的短距通道和短距和長距的穩(wěn)健版本給出了最佳結(jié)果浇冰。在穩(wěn)健模型中，僅使用長距通道仍然對測試統(tǒng)計的顯著性有大量過度報告聋亡。
最后肘习，圖9顯示了，使用修正的MVGC模型代替AR偏相關(guān)性進行了相同的分析坡倔。圖9(a)顯示了無運動模擬數(shù)據(jù)的結(jié)果漂佩。結(jié)果發(fā)現(xiàn)MVGC對模型中使用的通道不太敏感，短距罪塔、長距和甚至是兩種通道模型都具有相似的性能投蝉。據(jù)觀察，即使在靜止數(shù)據(jù)模擬的情況下征堪，MVGC的穩(wěn)健版本對測試統(tǒng)計數(shù)據(jù)的顯著性也有大量的過度報告瘩缆。研究者認為，這是模型中關(guān)于估計有效自由度的未修正誤差佃蚜。特別是庸娱，由于測量變量y着绊、x和/或Z中的時間點都可能是運動引起的異常值，因此難以在穩(wěn)健估計器中定義統(tǒng)計異常值熟尉。該模型中使用的Tukey雙平方權(quán)重假設(shè)離群值是獨立的归露，但當(dāng)前的模型沒有有效地考慮多個通道的共享離群值時間點，因此自由度稍微高估斤儿。在圖9(b)中剧包，針對添加運動偽影的數(shù)據(jù)模擬比較MVGC模型，結(jié)果發(fā)現(xiàn)穩(wěn)健方法優(yōu)于非穩(wěn)健方法往果。
建模滯后連接
修正后的MVGC同時考慮了零時滯和非零時滯疆液。因此，當(dāng)僅使用零時滯連接模擬數(shù)據(jù)時棚放，AR偏相關(guān)性被確定為最佳腾么。標(biāo)準(zhǔn)MVGC（非零時滯）在前幾節(jié)中沒有出現(xiàn)梅垄，是因為它在零非零時滯模擬中沒有發(fā)現(xiàn)任何結(jié)果[圖10(a)]。同樣艳吠，當(dāng)數(shù)據(jù)僅以非零時滯項的關(guān)系進行模擬時福也，AR偏相關(guān)模型失敗局骤，如圖10(b)所示。修正后的MVGC是唯一一個在兩種情況下都有效的模型暴凑。
在圖11中峦甩，在添加了運動偽影的情況下，顯示出了相同的模式现喳，其中ROC曲線對于零時滯模擬的情況下的穩(wěn)健AR偏相關(guān)和穩(wěn)健修正MVGC模型是最佳的凯傲，并且對于非零時滯模擬下的穩(wěn)健MVGC和穩(wěn)健修正的MVGC是最佳的。在這兩種情況下嗦篱，修正后的MVGC都不如AR偏相關(guān)或MVGC好冰单，因此如果知道存在哪種連接，則修正后的MVGC就不是首選灸促。反之诫欠，這是兩種情況下唯一可行的模型。
靜息態(tài)fNIRS數(shù)據(jù)的結(jié)果
為了檢驗本研究的功能連接方法對實驗靜息態(tài)數(shù)據(jù)的適用性浴栽，研究者比較了Pearson相關(guān)荒叼、AR相關(guān)和AR偏相關(guān)的穩(wěn)健版本，以在24名參與者樣本上模擬通道之間的功能連接典鸡。圖12顯示了三種方法連接路徑的組水平結(jié)果（q＜0.001被廓，F(xiàn)DR校正）。結(jié)果表明萝玷，與AR相關(guān)和Pearson相關(guān)相比嫁乘，使用AR偏相關(guān)導(dǎo)致更稀疏的網(wǎng)絡(luò)英遭，因為它校正了全局生理信號和自相關(guān)。使用Pearson相關(guān)性揭示了一個非常密集的連接網(wǎng)絡(luò)亦渗，但大多數(shù)連接可能是假陽性的挖诸。雖然沒有已知的地面真值或黃金標(biāo)準(zhǔn)的RS-fMRI網(wǎng)絡(luò)可用于定量比較這些方法，但AR偏相關(guān)可能是三種方法中假陽性率最低的一種方法法精。
結(jié)論
復(fù)雜的數(shù)值模擬模型在理論上可以獲取更多關(guān)于fNIRS時間通道之間關(guān)系的信息多律，但在實踐中很難求解。盡管修正的格蘭杰因果關(guān)系等方法可以對零時滯和時滯關(guān)系進行建模搂蜓，但當(dāng)連通性信息存在零時滯時狼荞，預(yù)白化偏相關(guān)性表現(xiàn)最佳，而當(dāng)時間序列之間存在時滯關(guān)系時帮碰，格蘭杰因果關(guān)系模型表現(xiàn)最佳相味。修正的MVGC在兩種情況下都表現(xiàn)良好，但在低信噪比情況下殉挽，自由度較低丰涉，靈敏度可能更低。同樣斯碌，帶通濾波雖然可能會提高信噪比和靈敏度一死，但也可能導(dǎo)致自由度損失和自相關(guān)增加，從而抵消增益傻唾。因此投慈，建議仔細考慮預(yù)處理流程中的步驟，以最大限度地提高靈敏度冠骄，同時將假陽性率降低到期望水平伪煤。
原文：Correction of global physiology in resting-state functional near-infrared spectroscopy.
DOI: 10.1117/1.NPh.9.3.035003
小伙伴們點個“在看”，加
(星標(biāo))關(guān)注茗創(chuàng)科技凛辣，將第一時間收到精彩內(nèi)容推送哦～
