傳統基于腦電（EEG）的抑郁癥研究將電極視為孤立節點，忽略了它們之間的關聯性，難以發掘抑郁癥患者異常大腦拓撲改變。為此，本文提出一種基于腦功能網絡（BFN）的抑郁癥識別框架，為避免容積導體效應，相位延遲指數用于構建BFN；以加權與二值化BFN信息互補為基礎，選取“小世界”特性密切相關及最小生成樹特定腦區BFN指標，采用遞進式指標分析策略尋找抑郁癥識別潛在標識物。本文以48名受試者靜息態EEG數據用于驗證方案，結果表明組間同步性在左顳、右頂枕、右額腦區明顯改變；加權BFN最短路徑長度和聚類系數，二值化BFN左顳和右額的葉子分數及右頂枕的直徑與患者健康問卷9項（PHQ-9）之間具有相關性，且獲得最高94.11%的識別率。此外，研究發現相對于健康對照者，抑郁癥患者的信息處理能力明顯下降。通過上述結論，期望本研究結果可為BFN構建與分析提供新的思路，為抑郁癥識別潛在標識物的發掘提供新的方法。
引用本文： 張冰濤, 周文穎, 李延林, 常文文, 徐斌斌. 基于腦功能網絡的抑郁癥識別研究. 生物醫學工程學雜志, 2022, 39(1): 47-55. doi: 10.7507/1001-5515.202108034 復制
引言
抑郁癥是一種常見精神疾病，持續情緒低落、興趣減退是其典型臨床癥狀[1]。據世界衛生組織披露，全球約有2.64億人不同程度受其影響[2]，每年約有85萬人直接或間接因其導致自殺身亡。重度抑郁癥復發率約為78%左右，完全治愈概率小[3]。相對于多數身體疾病，抑郁癥潛在神經機制和病理學機制尚不清楚。針對抑郁癥高發病率、高死亡率、高復發率等特點，為了避免隨時間推移造成大腦不可逆轉的損傷，本文基于腦功能網絡（brain functional network，BFN）展開了抑郁癥識別研究。
近年來，隨著生物醫學工程、神經科學、計算機科學等多學科交叉融合，研究人員發現大腦生理信號的異常改變與抑郁癥密切相關，他們以磁共振成像（magnetic resonance imaging，MRI）、腦磁圖（magnetoencephalography，MEG）和腦電圖（electroencephalograms，EEG）為基礎數據分別從大腦結構和功能視角展開了大量相關研究。Yao等[4]通過對抑郁癥患者MRI形態學特征研究發現，隨著抑郁癥程度不斷加深患者海馬區和杏仁核區的萎縮不斷向其他區域擴展。基于MEG特征Alageel等[5]研究發現抑郁癥患者右頂葉和額葉皮質活動明顯弱化。更多研究以EEG信號特征為切入點，先后利用多模態EEG信號特征融合、EEG信號特征與語音信號特征融合、EEG信號特征不同生成策略實現了自動抑郁癥識別[3, 6-7]。
基于特征的抑郁癥識別方法有其固有優勢，且人腦又是一種極為復雜的系統，若僅以形態學特征或生理學特征為基礎開展研究，難以揭示抑郁癥患者異常的大腦拓撲改變規律。隨著以圖論為基礎的復雜網絡發展，使用MRI、MEG和EEG等生理數據構建BFN成為可能，隨著對BFN核心拓撲研究的不斷深入，越來越多證據表明BFN拓撲異常改變與抑郁癥存在密切關聯[8]。Yu等[9]指出連通性是研究精神疾病患者腦功能異常的有效方法，BFN為神經系統疾病識別研究提供了新的途徑。
從BFN層面對抑郁癥識別進行研究，需解決兩個關鍵問題：① 何種類型數據更適合于BFN構建？② 如何分析BFN以發掘抑郁癥識別的潛在標志物。針對第一個關鍵問題，眾所周知，MRI數據具有較高的空間分辨率，MEG數據和EEG數據具有較高的時間分辨率[10]。BFN連通性被定義為空間上分離的不同神經元之間協調的時間依賴性活動[11]。因此，相對于MRI數據，EEG和MEG數據更適合于BFN構建。且EEG數據獲取具有非侵入、安全、易采集、低成本等優勢，經綜合考慮，本文使用EEG數據作為BFN構建的基礎數據。而第二個關鍵問題，本文則以加權與二值化BFN信息互補為基礎展開研究與探索。因此，本研究的開展將為BFN構建與分析提供新思路，為抑郁癥識別潛在標識物的發掘提供新方法。
1 實驗及數據預處理
1.1 實驗數據
本研究所用EEG數據來自精神障礙分析的多模式開放數據集（a multi-modal open dataset for mental-disorder analysis，MODMA Dataset）（http://modma.lzu.edu.cn/data/index/）[12]，本課題組已簽署實驗數據使用協議，并獲得授權可以使用該數據集中的全部數據。MODMA Dataset包含四種類型的抑郁癥患者及對應健康對照受試者數據，包括：① 128通道的事件相關電位（event-related potentials，ERP）數據；② 128通道的靜息態EEG數據；③ 3通道的靜息態EEG數據；④ 音頻數據。
本研究使用MODMA Dataset中第二類數據，該數據包括24名抑郁癥患者和29名健康對照受試者，基本信息如下：① 抑郁癥患者，患者健康問卷9項（patient health questionnaire 9-items，PHQ-9）得分≥5分，簡明國際神經精神障礙訪談檢查（mini international neuropsychiatric interview，MINI）得分≥17分；② 健康對照受試者，PHQ-9得分<5分，MINI得分<17分；③ 受試者年齡18～52歲之間，小學以上文化程度，智力正常。④ EEG數據采樣率為250 Hz，采集時長約5 min。
為保證組間性別及年齡盡可能匹配，本研究僅使用了24名抑郁癥患者（男性∶女性 = 13∶11）和24名健康對照受試者（男性∶女性 = 13∶11）的EEG數據。兩組PHQ-9得分的均值±標準差分別為：18.33 ± 3.50和2.58 ± 1.79，組間P < 0.001，差異具有統計學意義。在保證研究有效性和減少計算量的基礎上，本文選取均勻分布64通道EEG數據用于BFN的構建。
1.2 數據預處理
EEG數據采集過程不可避免地會引入噪聲，研究表明抑郁癥相關的EEG信號主要分布在0.5～50 Hz之間[13]，將有限沖激響應帶通濾波器用于消除低頻漂移和高頻噪聲。EEG信號噪聲主要來源于眼電（electrooculogram，EOG）干擾，即使靜息狀態下，眨眼所產生EOG信號振幅也是EEG信號振幅十倍以上，在0～16 Hz之間EOG信號噪聲與EEG信號嚴重重疊。為此，本文采用Peng等[14]提出的離散小波變換和自適應噪聲消除模型剔除EEG信號中的EOG信號噪聲，以獲得更加純凈的EEG數據。
2 本文方法
2.1 BFN構建
網絡G = (V, E)是由節點和節點之間的邊構成，V是節點集，E是邊集。本文中EEG信號對應電極被定義為節點，節點之間的依賴關系被定義為邊，同步性被作為依賴關系的載體。為了避免容積導體效應，相位延遲指數（phase lag index, PLI）用于計算節點之間的同步性[15]。任意EEG信號 
和 
，t時刻相位差 
計算過程如式（1）所示：
 |
其中，n和m為整數，神經科學中取值為1；
和 
是信號 
和 
的瞬時相位，如式（2）所示：
 |
其中，arctan是反正切函數，
是 
的希爾伯特變換，如式（3）所示：
 |
其中，PV是柯西主值。長度L的兩個時間序列之間PLI被定義為相位差 
的不對稱測量，如式（4）所示：
 |
其中，sign是符號函數。
研究表明抑郁癥患者與健康對照受試者的EEG信號在theta（4～8 Hz）頻段存在明顯差異[16]，EEG信號周期為12 s時，基于PLI構建的BFN趨于最佳狀態[17]。因此，本文使用theta頻段作為興趣頻段，每12 s計算一次EEG信號通道之間的PLI以生成鄰接矩陣。通過計算鄰接矩陣對應位置平均值得到64 × 64鄰接矩陣Cij，其定義如式（5）所示，進而繪制兩組的BFN。
 |
2.2 加權BFN指標
復雜網絡理論闡明了健康人群BFN具有穩定的“小世界”特性，即較小的最短路徑長度和較大的聚類系數，而此結論對于精神疾病患者并不成立[18]。因此，本文選擇最短路徑長度和聚類系數作為加權BFN的初始指標。
最短路徑長度能夠刻畫網絡全局連接特性和功能集成特性。若任意兩節點i和j之間存在最短路徑Lij，則加權BFN中任意兩個節點之間的最短路徑長度為Lij = 1/wij，其中wij是鄰接矩陣Cij中相應位置的元素，N是加權BFN的總節點數量。加權BFN的最短路徑長度L是節點最短路徑長度的平均值，定義如式（6）所示：
 |
聚類系數能夠刻畫網絡局部連接特性和功能隔離特性，節點聚類系數Ci是鄰接節點互為鄰居的概率，定義如式（7）所示：
 |
其中，i, k, j是加權BFN中任意三節點，wij, wik, wjk是式（5）鄰接矩陣中相應位置的元素。加權BFN的聚類系數C是節點聚類系數的平均值，定義如式（8）所示：
 |
2.3 二值化BFN指標
為研究組間差異性、提高組間可比性、避免加權BFN稀疏敏感度問題，本文以信息互補為理論基礎，采用最小生成樹（minimum spanning tree，MST）對加權BFN進行二值化。MST保證了組間差異性主要由拓撲差異產生，本文采用克魯斯卡爾算法構造MST，其過程如圖1所示，其中使用G = (V, E)表示原始BFN，V是節點集，E是邊集；MST(V, TE)是原始BFN對應的初始MST，TE是MST的邊集，初始值為空，即TE={}；ei是BFN中的任意邊。
圖1
加權BFN的MST構造
Figure1.
MST construction of weighted BFN
MST提供了原始網絡的核心拓撲信息[19]，本文選擇神經科學領域所公認的能夠有效分析精神疾病的典型MST屬性：度、葉子分數、直徑、介子中心、樹層次、離心率等作為二值化BFN的初始指標[20]。度k描述了與節點直接相連接的邊數，計算過程如式（9）所示：
 |
其中，v和u是節點集V中的任意兩個節點，evu表示由節點v和u所連接構成的邊。
葉子分數Lf描述了葉子節點數量與MST總節點數量之比，計算過程如式（10）所示：
 |
其中，L表示葉子節點數量，N表示MST總節點數量。
直徑d描述了給定節點的最短路徑長度，計算過程如式（11）所示：
 |
其中，TE(v)表示MST邊集TE中經過節點v邊的子集，max{}表示求集合中最大值。
介子中心BC描述了通過給定節點最短路徑的比例，計算過程如式（12）所示：
 |
其中，s, v, t是節點集V中的任意三個不相同節點，
是s和t之間最短路徑的數目，
是節點s和t之間最短路徑通過節點v的數目。
樹層次Th是MST的整合能力和中心節點過載之間平衡關系的量化度量，計算過程如式（13）所示：
 |
其中，BCmax是最大的BC。
離心率E描述了給定節點到其余節點距離的最大值，計算過程如式（14）所示：
 |
其中，d(v, u)表示v節點到節點集V中任意其它節點u的距離。
2.4 遞進式潛在標識物挖掘
為了發掘抑郁癥識別的潛在標識物，提出遞進式初始指標分析策略。第一步：采用非參數置換檢驗評估BFN組間差異具有統計學意義的初始指標[21]，首先計算BFN組間初始指標均值差異（t值）作為檢驗統計量的觀察值，接著全部受試者被隨機分配至抑郁癥組和健康對照組，保持組間受試者數量不變，重復計算1 000次重排列后組間t值，得到組間差異檢驗零分布，使用t值大于檢驗統計量觀察值的樣本排列比例作為組間差異P值，P < 0.05，則差異具有統計學意義。
基于PHQ-9的面對面訪談是抑郁癥診斷的基礎，生物指標或者BFN指標與量表得分之間相關性被廣泛用于評估相應指標是否能夠作為診斷有關疾病的潛在標識物[22-23]。基于此，遞進式分析策略第二步為：基于皮爾遜相關系數分析差異具有統計學意義的指標與PHQ-9得分之間的相關性，P < 0.05，則具有相關性。
2.5 識別性能評估
評估抑郁癥識別潛在標識物是否有效，本文采用K近鄰 (K-nearest neighbor，KNN)分類器，已有研究表明KNN分類器在抑郁癥識別領域優于其他分類器[6, 24-26]。為獲得具有統計學意義的結果，KNN分類器執行10折交叉驗證。本研究中準確率（accuracy），敏感性（sensitivity）和特異性（specificity）用于識別性能的評估。數學軟件MATLAB R2018a版（MathWorks Inc.，美國）用于實現KNN分類和10折交叉驗證。
3 結果
3.1 二值化BFN差異
克魯斯卡爾算法所生成兩組的鄰接矩陣及對應MST電極分布如圖2所示，此圖并不能完全清晰展示組間二值化BFN差異。為此，計算兩組鄰接矩陣的差異矩陣，進而繪制差異矩陣對應BFN的三維圖，如圖3所示，上排為左半球和右半球的外側視圖，下排為左半球和右半球的內側視圖，中間為全腦背側圖。從圖3左、右圖對應位置的三維圖可以看出，相比于健康對照組，抑郁癥組左顳（left temporal，LT）和右頂枕（right parietal-occipital，RPO）腦區的連接密度增加，即同步性增加；右額（right frontal，RF）腦區連接密度減少，即同步性減少。此結果其意味著抑郁癥組和健康對照組之間的二值化BFN差異主要分布在LT、RPO、RF三個腦區，即抑郁癥患者BFN在這三個腦區出現了拓撲異常。
圖2
鄰接矩陣及對應MST電極分布
Figure2.
Adjacency matrix and corresponding MST electrode distribution
圖3
差異矩陣對應BFN的三維圖
Figure3.
Three-dimensional diagram of difference matrix corresponding BFN
3.2 抑郁癥識別潛在標識物
本文以非參數置換檢驗對全腦加權BFN功能集成屬性最短路徑長度、功能隔離屬性聚類系數、以及二值化BFN的LT、RPO、RF等腦區初始指標進行評估。全腦加權BFN的最短路徑長度和聚類系數的組間差異結果如圖4所示，箱形圖中心線代表均值，星號（*）代表差異具有統計學意義（P < 0.05）。從圖4可看出，無論最短路徑長度還是聚類系數兩組間差異均具有統計學意義。相比于健康對照組，抑郁癥組最短路徑長度增加、聚類系數降低，意味著抑郁癥患者這兩項指標出現了異常，因此全腦加權BFN的最短路徑長度和聚類系數可作為抑郁癥識別的潛在標志物候選對象。二值化BFN在LT、RPO、RF三個腦區初始指標的組間差異結果如圖5所示，可以看出兩組間在LT腦區的葉子分數、樹層次、離心率，RPO腦區的直徑，以及RF腦區的葉子分數和直徑的差異均具有統計學意義，此意味這六項指標可作為抑郁癥識別潛在標志物的候選對象。
圖4
全腦加權BFN差異性分析（*P < 0.05）
Figure4.
Difference analysis of whole brain weighted BFN (*P < 0.05)
圖5
二值化BFN差異性分析（*P < 0.05）
Figure5.
Difference analysis of binarized BFN (*P < 0.05)
基于遞進式初始指標分析策略，使用皮爾遜相關系數評估具有統計學意義的BFN指標與PHQ-9得分之間相關性，結果如圖6所示，r表示相關系數，P值為概率，斜線表示擬合曲線。從圖6可知，全腦加權BFN最短路徑長度和聚類系數與PHQ-9得分兩者具有相關性，二值化BFN中LT和RF的葉子分數以及RPO的直徑與PHQ-9得分之間也具有相關性。因此，這五項BFN指標被作為抑郁癥識別的潛在標志物。考慮篇幅限制，圖中未列出與PHQ-9得分不相關的指標。
圖6
PHQ-9評分相關的BFN指標
Figure6.
BFN index correlated with PHQ-9 score
3.3 潛在標志物識別性能
評估上述五項BFN指標識別抑郁癥的性能，將其輸入至KNN分類器（k=3, 搜索算法=“kd-tree”）進行10折交叉驗證，當五項指標共同作為輸入時，最低分類準確率是86.42%，敏感性、特異性分別是81.56%和91.58%（如圖7所示：迭代-5）；最高分類準確率是94.35%，敏感性、特異性分別是91.84%和97.36%（如圖7所示：迭代-9）；10折交叉驗證的平均準確率是89.91%。因此，這五項指標共同作為潛在標識物可以有效識別抑郁癥。
圖7
10折交叉驗證的識別結果
Figure7.
Recognition result of 10-fold cross validation
每項指標單獨作為潛在標志物的抑郁癥識別性能結果如表1所示，五項指標作為潛在標志物識別能力由強到弱依次為：全腦聚類系數、最短路徑長度、LT-葉子分數、RF-葉子分數、RPO-直徑。
表1
BFN指標單獨作為潛在標識物的識別能力
Table1.
Recognition ability of BFN index alone as potential bio marker
3.4 泛化性能
為了驗證本文提出方法的泛化能力，采用來自亞利桑那大學公開抑郁癥EEG數據集（https://unmm-my.sharepoint.com/personal/jcavanagh_unm_edu/_layouts/15/onedrive.aspx）[27]進行驗證。該數據集共包括46名（男性∶女性 = 12∶34）抑郁癥患者和75名（男性∶女性 = 35∶40）健康對照受試者的64通道EEG數據。兩組的貝克抑郁自評量表（Beck depression inventory，BDI）均值 ± 標準差分別為22.22 ± 4.90和1.73 ± 1.65，組間P<0.01，差異具有統計學意義。
根據本文所提出方法，計算上述公開數據集中全腦加權BFN最短路徑長度和聚類系數、LT-葉子分數、RPO-直徑、RF-葉子分數，然后將這些指標共同輸入至KNN分類器（參數設置與上文相同）進行10折交叉驗證，獲得了87.53%的平均準確率，此結果略低于MODMA Dataset 上的89.91%平均準確率。此外，將五項指標分別作為抑郁癥識別的潛在標識物，識別結果如表2所示。可以看出無論五項潛在標識物共同或者獨立作為KNN分類器的輸入，準確率均有略微下降，但整體而言沒有明顯下降，說明了本文所提出方法在一定程度上具有健壯的泛化能力。
表2
五項指標泛化性能驗證結果
Table2.
Generalization performance verification result of five index
4 討論
由于研究方法的不同，有關抑郁癥的研究也存在相悖結論。Zhou等[28]發現抑郁癥患者大腦功能連接數量增加，而Connolly等[29]發現抑郁癥患者大腦功能連接數量減少。本文基于MST差異矩陣分析發現相對于健康對照受試者，抑郁癥患者BFN的同步性發生改變，LT和RPO腦區同步性明顯增加，RF腦區同步性明顯降低。作者推測此結果是抑郁癥患者大腦信息處理功能異常在BFN層面的表現形式，此結論也得到了相關研究的間接驗證，Joormann等[30]指出抑郁癥患者情緒調節能力存在明顯下降的趨勢，Leistedt等[31]發現抑郁癥患者BFN信息處理能力下降。本研究揭示抑郁癥患者部分腦區的BFN同步性存在失衡，但是抑郁癥患者大腦信息處理功能異常的原因還需進一步探索以獲得更為可靠的結論。
小世界網絡是一種最優的網絡組織結構，健康人群EEG數據或MEG數據構建的BFN具有穩定的“小世界”特征，而通常精神疾病患者的BFN“小世界”特征往往出現了弱化趨勢[18]。本研究發現抑郁癥組與健康對照組的最短路徑長度和聚類系數差異具有統計學意義，全腦加權BFN最短路徑長度增加，聚類系數減小。其意味著抑郁癥患者“小世界”特性存在隨機化趨勢，為進一步驗證該結論，本文分析了兩組的 “小世界”特征指數，其定義如式（15）所示：
 |
Creal和Lreal是原始網絡的聚類系數與最短路徑長度，Crandom和Lrandom是原始網絡對應隨機網絡的聚類系數與最短路徑長度。
> 1代表網絡具有“小世界”特性，且 
值越大網絡的“小世界”特性越強[32]。抑郁癥組和健康對照組的“小世界”特性指數差異性結果如圖8所示，可以看出兩組的
值均大于1，其說明了兩組的BFN均具有“小世界”特性。相比于健康對照組，抑郁癥組“小世界”特性指數明顯下降，該結果意味著抑郁癥患者BFN出現了隨機化趨勢。
圖8
組間“小世界”特征指數對比
Figure8.
Comparison of “small world” characteristic indexes between groups
5 結論
本文旨在通過探索抑郁癥患者BFN異常變化，進而不斷向抑郁癥發病神經機制真相逼近。基于靜息態EEG數據的加權和二值化BFN組間對比發現，相對于健康對照組，抑郁癥組BFN同步性發生了改變，如LT和RPO腦區同步性增加，RF腦區同步性減少，這些改變可能是抑郁癥患者大腦信息處理功能異常在BFN層面的外在表現形式。通過對比組間“小世界”特征指數發現抑郁癥患者的BFN出現了隨機化趨勢。全腦加權BFN的聚類系數作為KNN分類器輸入可獲得最高94.11%的抑郁癥識別率。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：張冰濤負責研究方案設計，文章構思，論文撰寫；周文穎負責實驗數據處理與論文修改；李延林負責數據分析；常文文負責論文指導，審閱；徐斌斌負責圖表繪制。
倫理聲明：本研究通過了蘭州大學第二附屬醫院倫理委員會的審批（批文編號：2015-017）。
引言
抑郁癥是一種常見精神疾病，持續情緒低落、興趣減退是其典型臨床癥狀[1]。據世界衛生組織披露，全球約有2.64億人不同程度受其影響[2]，每年約有85萬人直接或間接因其導致自殺身亡。重度抑郁癥復發率約為78%左右，完全治愈概率小[3]。相對于多數身體疾病，抑郁癥潛在神經機制和病理學機制尚不清楚。針對抑郁癥高發病率、高死亡率、高復發率等特點，為了避免隨時間推移造成大腦不可逆轉的損傷，本文基于腦功能網絡（brain functional network，BFN）展開了抑郁癥識別研究。
近年來，隨著生物醫學工程、神經科學、計算機科學等多學科交叉融合，研究人員發現大腦生理信號的異常改變與抑郁癥密切相關，他們以磁共振成像（magnetic resonance imaging，MRI）、腦磁圖（magnetoencephalography，MEG）和腦電圖（electroencephalograms，EEG）為基礎數據分別從大腦結構和功能視角展開了大量相關研究。Yao等[4]通過對抑郁癥患者MRI形態學特征研究發現，隨著抑郁癥程度不斷加深患者海馬區和杏仁核區的萎縮不斷向其他區域擴展。基于MEG特征Alageel等[5]研究發現抑郁癥患者右頂葉和額葉皮質活動明顯弱化。更多研究以EEG信號特征為切入點，先后利用多模態EEG信號特征融合、EEG信號特征與語音信號特征融合、EEG信號特征不同生成策略實現了自動抑郁癥識別[3, 6-7]。
基于特征的抑郁癥識別方法有其固有優勢，且人腦又是一種極為復雜的系統，若僅以形態學特征或生理學特征為基礎開展研究，難以揭示抑郁癥患者異常的大腦拓撲改變規律。隨著以圖論為基礎的復雜網絡發展，使用MRI、MEG和EEG等生理數據構建BFN成為可能，隨著對BFN核心拓撲研究的不斷深入，越來越多證據表明BFN拓撲異常改變與抑郁癥存在密切關聯[8]。Yu等[9]指出連通性是研究精神疾病患者腦功能異常的有效方法，BFN為神經系統疾病識別研究提供了新的途徑。
從BFN層面對抑郁癥識別進行研究，需解決兩個關鍵問題：① 何種類型數據更適合于BFN構建？② 如何分析BFN以發掘抑郁癥識別的潛在標志物。針對第一個關鍵問題，眾所周知，MRI數據具有較高的空間分辨率，MEG數據和EEG數據具有較高的時間分辨率[10]。BFN連通性被定義為空間上分離的不同神經元之間協調的時間依賴性活動[11]。因此，相對于MRI數據，EEG和MEG數據更適合于BFN構建。且EEG數據獲取具有非侵入、安全、易采集、低成本等優勢，經綜合考慮，本文使用EEG數據作為BFN構建的基礎數據。而第二個關鍵問題，本文則以加權與二值化BFN信息互補為基礎展開研究與探索。因此，本研究的開展將為BFN構建與分析提供新思路，為抑郁癥識別潛在標識物的發掘提供新方法。
1 實驗及數據預處理
1.1 實驗數據
本研究所用EEG數據來自精神障礙分析的多模式開放數據集（a multi-modal open dataset for mental-disorder analysis，MODMA Dataset）（http://modma.lzu.edu.cn/data/index/）[12]，本課題組已簽署實驗數據使用協議，并獲得授權可以使用該數據集中的全部數據。MODMA Dataset包含四種類型的抑郁癥患者及對應健康對照受試者數據，包括：① 128通道的事件相關電位（event-related potentials，ERP）數據；② 128通道的靜息態EEG數據；③ 3通道的靜息態EEG數據；④ 音頻數據。
本研究使用MODMA Dataset中第二類數據，該數據包括24名抑郁癥患者和29名健康對照受試者，基本信息如下：① 抑郁癥患者，患者健康問卷9項（patient health questionnaire 9-items，PHQ-9）得分≥5分，簡明國際神經精神障礙訪談檢查（mini international neuropsychiatric interview，MINI）得分≥17分；② 健康對照受試者，PHQ-9得分<5分，MINI得分<17分；③ 受試者年齡18～52歲之間，小學以上文化程度，智力正常。④ EEG數據采樣率為250 Hz，采集時長約5 min。
為保證組間性別及年齡盡可能匹配，本研究僅使用了24名抑郁癥患者（男性∶女性 = 13∶11）和24名健康對照受試者（男性∶女性 = 13∶11）的EEG數據。兩組PHQ-9得分的均值±標準差分別為：18.33 ± 3.50和2.58 ± 1.79，組間P < 0.001，差異具有統計學意義。在保證研究有效性和減少計算量的基礎上，本文選取均勻分布64通道EEG數據用于BFN的構建。
1.2 數據預處理
EEG數據采集過程不可避免地會引入噪聲，研究表明抑郁癥相關的EEG信號主要分布在0.5～50 Hz之間[13]，將有限沖激響應帶通濾波器用于消除低頻漂移和高頻噪聲。EEG信號噪聲主要來源于眼電（electrooculogram，EOG）干擾，即使靜息狀態下，眨眼所產生EOG信號振幅也是EEG信號振幅十倍以上，在0～16 Hz之間EOG信號噪聲與EEG信號嚴重重疊。為此，本文采用Peng等[14]提出的離散小波變換和自適應噪聲消除模型剔除EEG信號中的EOG信號噪聲，以獲得更加純凈的EEG數據。
2 本文方法
2.1 BFN構建
網絡G = (V, E)是由節點和節點之間的邊構成，V是節點集，E是邊集。本文中EEG信號對應電極被定義為節點，節點之間的依賴關系被定義為邊，同步性被作為依賴關系的載體。為了避免容積導體效應，相位延遲指數（phase lag index, PLI）用于計算節點之間的同步性[15]。任意EEG信號 和 ，t時刻相位差 計算過程如式（1）所示：
|
其中，n和m為整數，神經科學中取值為1； 和 是信號 和 的瞬時相位，如式（2）所示：
|
其中，arctan是反正切函數， 是 的希爾伯特變換，如式（3）所示：
|
其中，PV是柯西主值。長度L的兩個時間序列之間PLI被定義為相位差 的不對稱測量，如式（4）所示：
|
其中，sign是符號函數。
研究表明抑郁癥患者與健康對照受試者的EEG信號在theta（4～8 Hz）頻段存在明顯差異[16]，EEG信號周期為12 s時，基于PLI構建的BFN趨于最佳狀態[17]。因此，本文使用theta頻段作為興趣頻段，每12 s計算一次EEG信號通道之間的PLI以生成鄰接矩陣。通過計算鄰接矩陣對應位置平均值得到64 × 64鄰接矩陣Cij，其定義如式（5）所示，進而繪制兩組的BFN。
|
2.2 加權BFN指標
復雜網絡理論闡明了健康人群BFN具有穩定的“小世界”特性，即較小的最短路徑長度和較大的聚類系數，而此結論對于精神疾病患者并不成立[18]。因此，本文選擇最短路徑長度和聚類系數作為加權BFN的初始指標。
最短路徑長度能夠刻畫網絡全局連接特性和功能集成特性。若任意兩節點i和j之間存在最短路徑Lij，則加權BFN中任意兩個節點之間的最短路徑長度為Lij = 1/wij，其中wij是鄰接矩陣Cij中相應位置的元素，N是加權BFN的總節點數量。加權BFN的最短路徑長度L是節點最短路徑長度的平均值，定義如式（6）所示：
|
聚類系數能夠刻畫網絡局部連接特性和功能隔離特性，節點聚類系數Ci是鄰接節點互為鄰居的概率，定義如式（7）所示：
|
其中，i, k, j是加權BFN中任意三節點，wij, wik, wjk是式（5）鄰接矩陣中相應位置的元素。加權BFN的聚類系數C是節點聚類系數的平均值，定義如式（8）所示：
|
2.3 二值化BFN指標
為研究組間差異性、提高組間可比性、避免加權BFN稀疏敏感度問題，本文以信息互補為理論基礎，采用最小生成樹（minimum spanning tree，MST）對加權BFN進行二值化。MST保證了組間差異性主要由拓撲差異產生，本文采用克魯斯卡爾算法構造MST，其過程如圖1所示，其中使用G = (V, E)表示原始BFN，V是節點集，E是邊集；MST(V, TE)是原始BFN對應的初始MST，TE是MST的邊集，初始值為空，即TE={}；ei是BFN中的任意邊。
圖1
加權BFN的MST構造
Figure1.
MST construction of weighted BFN
MST提供了原始網絡的核心拓撲信息[19]，本文選擇神經科學領域所公認的能夠有效分析精神疾病的典型MST屬性：度、葉子分數、直徑、介子中心、樹層次、離心率等作為二值化BFN的初始指標[20]。度k描述了與節點直接相連接的邊數，計算過程如式（9）所示：
|
其中，v和u是節點集V中的任意兩個節點，evu表示由節點v和u所連接構成的邊。
葉子分數Lf描述了葉子節點數量與MST總節點數量之比，計算過程如式（10）所示：
|
其中，L表示葉子節點數量，N表示MST總節點數量。
直徑d描述了給定節點的最短路徑長度，計算過程如式（11）所示：
|
其中，TE(v)表示MST邊集TE中經過節點v邊的子集，max{}表示求集合中最大值。
介子中心BC描述了通過給定節點最短路徑的比例，計算過程如式（12）所示：
|
其中，s, v, t是節點集V中的任意三個不相同節點，是s和t之間最短路徑的數目，是節點s和t之間最短路徑通過節點v的數目。
樹層次Th是MST的整合能力和中心節點過載之間平衡關系的量化度量，計算過程如式（13）所示：
|
其中，BCmax是最大的BC。
離心率E描述了給定節點到其余節點距離的最大值，計算過程如式（14）所示：
|
其中，d(v, u)表示v節點到節點集V中任意其它節點u的距離。
2.4 遞進式潛在標識物挖掘
為了發掘抑郁癥識別的潛在標識物，提出遞進式初始指標分析策略。第一步：采用非參數置換檢驗評估BFN組間差異具有統計學意義的初始指標[21]，首先計算BFN組間初始指標均值差異（t值）作為檢驗統計量的觀察值，接著全部受試者被隨機分配至抑郁癥組和健康對照組，保持組間受試者數量不變，重復計算1 000次重排列后組間t值，得到組間差異檢驗零分布，使用t值大于檢驗統計量觀察值的樣本排列比例作為組間差異P值，P < 0.05，則差異具有統計學意義。
基于PHQ-9的面對面訪談是抑郁癥診斷的基礎，生物指標或者BFN指標與量表得分之間相關性被廣泛用于評估相應指標是否能夠作為診斷有關疾病的潛在標識物[22-23]。基于此，遞進式分析策略第二步為：基于皮爾遜相關系數分析差異具有統計學意義的指標與PHQ-9得分之間的相關性，P < 0.05，則具有相關性。
2.5 識別性能評估
評估抑郁癥識別潛在標識物是否有效，本文采用K近鄰 (K-nearest neighbor，KNN)分類器，已有研究表明KNN分類器在抑郁癥識別領域優于其他分類器[6, 24-26]。為獲得具有統計學意義的結果，KNN分類器執行10折交叉驗證。本研究中準確率（accuracy），敏感性（sensitivity）和特異性（specificity）用于識別性能的評估。數學軟件MATLAB R2018a版（MathWorks Inc.，美國）用于實現KNN分類和10折交叉驗證。
3 結果
3.1 二值化BFN差異
克魯斯卡爾算法所生成兩組的鄰接矩陣及對應MST電極分布如圖2所示，此圖并不能完全清晰展示組間二值化BFN差異。為此，計算兩組鄰接矩陣的差異矩陣，進而繪制差異矩陣對應BFN的三維圖，如圖3所示，上排為左半球和右半球的外側視圖，下排為左半球和右半球的內側視圖，中間為全腦背側圖。從圖3左、右圖對應位置的三維圖可以看出，相比于健康對照組，抑郁癥組左顳（left temporal，LT）和右頂枕（right parietal-occipital，RPO）腦區的連接密度增加，即同步性增加；右額（right frontal，RF）腦區連接密度減少，即同步性減少。此結果其意味著抑郁癥組和健康對照組之間的二值化BFN差異主要分布在LT、RPO、RF三個腦區，即抑郁癥患者BFN在這三個腦區出現了拓撲異常。
圖2
鄰接矩陣及對應MST電極分布
Figure2.
Adjacency matrix and corresponding MST electrode distribution
圖3
差異矩陣對應BFN的三維圖
Figure3.
Three-dimensional diagram of difference matrix corresponding BFN
3.2 抑郁癥識別潛在標識物
本文以非參數置換檢驗對全腦加權BFN功能集成屬性最短路徑長度、功能隔離屬性聚類系數、以及二值化BFN的LT、RPO、RF等腦區初始指標進行評估。全腦加權BFN的最短路徑長度和聚類系數的組間差異結果如圖4所示，箱形圖中心線代表均值，星號（*）代表差異具有統計學意義（P < 0.05）。從圖4可看出，無論最短路徑長度還是聚類系數兩組間差異均具有統計學意義。相比于健康對照組，抑郁癥組最短路徑長度增加、聚類系數降低，意味著抑郁癥患者這兩項指標出現了異常，因此全腦加權BFN的最短路徑長度和聚類系數可作為抑郁癥識別的潛在標志物候選對象。二值化BFN在LT、RPO、RF三個腦區初始指標的組間差異結果如圖5所示，可以看出兩組間在LT腦區的葉子分數、樹層次、離心率，RPO腦區的直徑，以及RF腦區的葉子分數和直徑的差異均具有統計學意義，此意味這六項指標可作為抑郁癥識別潛在標志物的候選對象。
圖4
全腦加權BFN差異性分析（*P < 0.05）
Figure4.
Difference analysis of whole brain weighted BFN (*P < 0.05)
圖5
二值化BFN差異性分析（*P < 0.05）
Figure5.
Difference analysis of binarized BFN (*P < 0.05)
基于遞進式初始指標分析策略，使用皮爾遜相關系數評估具有統計學意義的BFN指標與PHQ-9得分之間相關性，結果如圖6所示，r表示相關系數，P值為概率，斜線表示擬合曲線。從圖6可知，全腦加權BFN最短路徑長度和聚類系數與PHQ-9得分兩者具有相關性，二值化BFN中LT和RF的葉子分數以及RPO的直徑與PHQ-9得分之間也具有相關性。因此，這五項BFN指標被作為抑郁癥識別的潛在標志物。考慮篇幅限制，圖中未列出與PHQ-9得分不相關的指標。
圖6
PHQ-9評分相關的BFN指標
Figure6.
BFN index correlated with PHQ-9 score
3.3 潛在標志物識別性能
評估上述五項BFN指標識別抑郁癥的性能，將其輸入至KNN分類器（k=3, 搜索算法=“kd-tree”）進行10折交叉驗證，當五項指標共同作為輸入時，最低分類準確率是86.42%，敏感性、特異性分別是81.56%和91.58%（如圖7所示：迭代-5）；最高分類準確率是94.35%，敏感性、特異性分別是91.84%和97.36%（如圖7所示：迭代-9）；10折交叉驗證的平均準確率是89.91%。因此，這五項指標共同作為潛在標識物可以有效識別抑郁癥。
圖7
10折交叉驗證的識別結果
Figure7.
Recognition result of 10-fold cross validation
每項指標單獨作為潛在標志物的抑郁癥識別性能結果如表1所示，五項指標作為潛在標志物識別能力由強到弱依次為：全腦聚類系數、最短路徑長度、LT-葉子分數、RF-葉子分數、RPO-直徑。
表1
BFN指標單獨作為潛在標識物的識別能力
Table1.
Recognition ability of BFN index alone as potential bio marker
3.4 泛化性能
為了驗證本文提出方法的泛化能力，采用來自亞利桑那大學公開抑郁癥EEG數據集（https://unmm-my.sharepoint.com/personal/jcavanagh_unm_edu/_layouts/15/onedrive.aspx）[27]進行驗證。該數據集共包括46名（男性∶女性 = 12∶34）抑郁癥患者和75名（男性∶女性 = 35∶40）健康對照受試者的64通道EEG數據。兩組的貝克抑郁自評量表（Beck depression inventory，BDI）均值 ± 標準差分別為22.22 ± 4.90和1.73 ± 1.65，組間P<0.01，差異具有統計學意義。
根據本文所提出方法，計算上述公開數據集中全腦加權BFN最短路徑長度和聚類系數、LT-葉子分數、RPO-直徑、RF-葉子分數，然后將這些指標共同輸入至KNN分類器（參數設置與上文相同）進行10折交叉驗證，獲得了87.53%的平均準確率，此結果略低于MODMA Dataset 上的89.91%平均準確率。此外，將五項指標分別作為抑郁癥識別的潛在標識物，識別結果如表2所示。可以看出無論五項潛在標識物共同或者獨立作為KNN分類器的輸入，準確率均有略微下降，但整體而言沒有明顯下降，說明了本文所提出方法在一定程度上具有健壯的泛化能力。
表2
五項指標泛化性能驗證結果
Table2.
Generalization performance verification result of five index
4 討論
由于研究方法的不同，有關抑郁癥的研究也存在相悖結論。Zhou等[28]發現抑郁癥患者大腦功能連接數量增加，而Connolly等[29]發現抑郁癥患者大腦功能連接數量減少。本文基于MST差異矩陣分析發現相對于健康對照受試者，抑郁癥患者BFN的同步性發生改變，LT和RPO腦區同步性明顯增加，RF腦區同步性明顯降低。作者推測此結果是抑郁癥患者大腦信息處理功能異常在BFN層面的表現形式，此結論也得到了相關研究的間接驗證，Joormann等[30]指出抑郁癥患者情緒調節能力存在明顯下降的趨勢，Leistedt等[31]發現抑郁癥患者BFN信息處理能力下降。本研究揭示抑郁癥患者部分腦區的BFN同步性存在失衡，但是抑郁癥患者大腦信息處理功能異常的原因還需進一步探索以獲得更為可靠的結論。
小世界網絡是一種最優的網絡組織結構，健康人群EEG數據或MEG數據構建的BFN具有穩定的“小世界”特征，而通常精神疾病患者的BFN“小世界”特征往往出現了弱化趨勢[18]。本研究發現抑郁癥組與健康對照組的最短路徑長度和聚類系數差異具有統計學意義，全腦加權BFN最短路徑長度增加，聚類系數減小。其意味著抑郁癥患者“小世界”特性存在隨機化趨勢，為進一步驗證該結論，本文分析了兩組的 “小世界”特征指數，其定義如式（15）所示：
|
Creal和Lreal是原始網絡的聚類系數與最短路徑長度，Crandom和Lrandom是原始網絡對應隨機網絡的聚類系數與最短路徑長度。 > 1代表網絡具有“小世界”特性，且 值越大網絡的“小世界”特性越強[32]。抑郁癥組和健康對照組的“小世界”特性指數差異性結果如圖8所示，可以看出兩組的值均大于1，其說明了兩組的BFN均具有“小世界”特性。相比于健康對照組，抑郁癥組“小世界”特性指數明顯下降，該結果意味著抑郁癥患者BFN出現了隨機化趨勢。
圖8
組間“小世界”特征指數對比
Figure8.
Comparison of “small world” characteristic indexes between groups
5 結論
本文旨在通過探索抑郁癥患者BFN異常變化，進而不斷向抑郁癥發病神經機制真相逼近。基于靜息態EEG數據的加權和二值化BFN組間對比發現，相對于健康對照組，抑郁癥組BFN同步性發生了改變，如LT和RPO腦區同步性增加，RF腦區同步性減少，這些改變可能是抑郁癥患者大腦信息處理功能異常在BFN層面的外在表現形式。通過對比組間“小世界”特征指數發現抑郁癥患者的BFN出現了隨機化趨勢。全腦加權BFN的聚類系數作為KNN分類器輸入可獲得最高94.11%的抑郁癥識別率。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：張冰濤負責研究方案設計，文章構思，論文撰寫；周文穎負責實驗數據處理與論文修改；李延林負責數據分析；常文文負責論文指導，審閱；徐斌斌負責圖表繪制。
倫理聲明：本研究通過了蘭州大學第二附屬醫院倫理委員會的審批（批文編號：2015-017）。
