柯　銘　曹　黎

(蘭州理工大學計算機與通信學院　甘肅 蘭州 730050)

0　引　言

復雜網(wǎng)絡理論以Watts和Strogatz關于小世界網(wǎng)絡[1]，以及Barabasi和Albert關于無標度網(wǎng)絡[2]的開創(chuàng)性工作之后，已經(jīng)發(fā)展成為一個多學科交叉研究的前沿。復雜網(wǎng)絡結合腦的研究已經(jīng)成為腦科學領域研究的熱點。最早，Stam[3]利用5個健康被試在無任務條件下的腦磁圖數(shù)據(jù)，建立了人腦功能網(wǎng)絡并研究了不同頻率段下構建的腦網(wǎng)絡的網(wǎng)絡屬性。隨后，Salvador等[4-6]進行了一系列基于fMRI數(shù)據(jù)構造的復雜網(wǎng)絡分析，其研究的主要結果為：進行了靜息態(tài)fMRI的頻率域分析，發(fā)現(xiàn)兩兩腦區(qū)之間的部分相干譜取決于腦區(qū)間的解剖距離；發(fā)現(xiàn)靜息態(tài)fMRI的低頻條件下的腦網(wǎng)絡具有“小世界”屬性；構建了老年與青年的靜息態(tài)fMRI的腦網(wǎng)絡，發(fā)現(xiàn)兩個網(wǎng)絡均呈現(xiàn)“小世界”屬性以及較高的全局效率和局部效率等結果。Basset等[7]綜述了近些年該領域的發(fā)現(xiàn)并指出“小世界”是適應性復雜動力學問題的一個解決方案。同時，王艷群等[8]構建了正常人的fMRI腦功能網(wǎng)絡，驗證了其具有模塊化結構。復雜網(wǎng)絡結合腦的研究同樣被應用到了疾病的研究。李越等[9]分析了抑郁癥患者腦網(wǎng)絡社團結構的改變，為該疾病提供了重要的應用與臨床診斷價值。Pedersen等[10]發(fā)現(xiàn)癲癇患者腦網(wǎng)絡的節(jié)點以更離散化的方式變化可能與該疾病相關。

目前，基于fMRI構建的復雜腦網(wǎng)絡旨在研究初始網(wǎng)絡的靜態(tài)拓撲屬性等指標。如上文所述的研究成果等，學者們較多關注的是節(jié)點并未遭到攻擊時的理想狀態(tài)下，腦網(wǎng)絡的各項網(wǎng)絡拓撲屬性。但腦網(wǎng)絡是一個高效的動態(tài)的網(wǎng)絡[7]，某個節(jié)點的失效可能會引發(fā)災難性的后果。尤其負荷最大節(jié)點的失效，其承載的信息量將瞬間崩潰，對整個腦網(wǎng)絡的魯棒性將會產(chǎn)生的影響尚不清楚。從保護腦網(wǎng)絡整體性角度，我們認為研究腦網(wǎng)絡負荷最大節(jié)點失效引起的動態(tài)變化極其重要。因此，本文引入一種相繼故障模型來研究網(wǎng)絡中負荷最大節(jié)點失效所引起整個網(wǎng)絡的動態(tài)變化。

相繼故障是指在大部分實際網(wǎng)絡中，一個或多個節(jié)點的失效(隨機或蓄意)會通過節(jié)點之間的耦合關系引起其他節(jié)點發(fā)生故障，由此產(chǎn)生連鎖反應，導致相當一部分節(jié)點甚至整個網(wǎng)絡的崩潰。CLM相繼故障模型[11]為負荷-容量模型中的一種子模型，亦稱為節(jié)點與邊混合動態(tài)模型。在該模型下，節(jié)點過載不是簡單地移除，而是降低該節(jié)點與其他節(jié)點之間的連接效率，與腦網(wǎng)絡較為符合。Kinney等[12]曾將CLM相繼故障模型應用到北美電力網(wǎng)絡。發(fā)現(xiàn)攻擊負荷最大的發(fā)電節(jié)點或者輸電節(jié)點時，網(wǎng)絡的整體性能相對于正常狀態(tài)下有所下降，且隨著容量系數(shù)的減小，網(wǎng)絡的平均效率亦下降。

本文首次將CLM相繼故障模型應用到正常人腦進行研究。我們旨在關注正常人腦網(wǎng)絡在節(jié)點失效后網(wǎng)絡拓撲結構的變化，探索網(wǎng)絡效率的改變，以此研究正常人腦網(wǎng)絡的魯棒性與脆弱性，具有創(chuàng)新意義。首先對正常人的靜息態(tài)功能核磁共振圖像進行預處理，利用AAL模板將正常人的大腦分割成90個腦區(qū)。然后分別計算兩兩腦區(qū)之間的Pearson相關系數(shù)。接著設定閾值以獲得鄰接矩陣來描述腦功能網(wǎng)絡。最后設置容量系數(shù)，賦予該網(wǎng)絡每個腦區(qū)容許值。并對負荷最大腦區(qū)進行攻擊，分析人腦功能網(wǎng)絡的魯棒性與脆弱性。

1　方　法

1.1　被　試

采用18例正常人的核磁共振(fMRI)圖像作為被試。其中，7例男性，11例女性，年齡在10～34歲之間，平均年齡19.9歲，且均為右利手。該數(shù)據(jù)是在Siemens Magnetom Trio 3.0T超導MR儀上進行靜息態(tài)功能磁共振數(shù)據(jù)采集。采用梯度回波-平面回波成像GRE-EPI(gradient echo-echo planar imaging)序列進行掃描，掃描范圍包括全腦，基線平行于前后聯(lián)合，掃描參數(shù)如下：TR 2 000 ms，TE 30 ms，層厚3.8 mm，層間距0.38 mm，層數(shù)34，F(xiàn)OV 240×240 mm ，矩陣64×64，反轉角90°，共采集200個時間點，掃描時間6 min 40 s 。要求受試者平臥，頭部固定，閉目，塞耳，盡量不做意向性思維。

1.2　數(shù)據(jù)預處理

在matlab環(huán)境下利用DPARSF和rest軟件對數(shù)據(jù)進行處理。其中，DPARSF進行的處理步驟包括：DICOM數(shù)據(jù)的格式轉換、去除前10個時間點、頭動校正(平動1 mm，旋轉1度)、空間標準化、平滑、去線性漂移、低頻濾波。利用rest提取協(xié)變量、再去除協(xié)變量。其中，協(xié)變量定義的是灰質(zhì)、白質(zhì)和腦脊液。

1.3　功能連接

腦區(qū)間的時間序列同步性可表示腦區(qū)間的功能關系。功能連接關系主要由Pearson相關系數(shù)[13]表示。本文Pearson相關系數(shù)旨在計算兩兩腦區(qū)間的時間相關性，以此表達功能連接關系。由此，采用rest軟件做功能連接。首先，利用AAL模板(共90個腦區(qū))定義ROI，將大腦分割成90個節(jié)點，并提取90個節(jié)點的時間序列。然后，利用rest軟件的“ROI-wise”功能，計算90個腦區(qū)間的Pearson相關系數(shù)，可得一90維對稱矩陣，值均在-1～1之間且包括端點值。該矩陣即表示所有腦區(qū)之間的Pearson相關系數(shù)。Pearson相關系數(shù)的計算公式為：

(1)

1.4　腦網(wǎng)絡構建

由功能連接，得到每個被試90個腦區(qū)兩兩之間的Pearson相關系數(shù)矩陣。首先，Z值化每個被試該相關系數(shù)矩陣，使其滿足正態(tài)性。然后，將所有被試的相關系數(shù)矩陣求和并均值化，得到一個表示所有被試的相關系數(shù)矩陣。最后，構建二值網(wǎng)絡時，為消除由噪聲引起的弱連接，需設置閾值a(a>0)將連接強度偏弱的值刪除。若相關系數(shù)矩陣中的值大于等于a或小于等于-a，則表示該兩腦區(qū)存在連接，即設置為1，否則設為0。鄰接矩陣中對角線的元素設置為0，以避免自連。據(jù)文獻[6]報導，腦網(wǎng)絡是一個復雜且稀疏的以大量短程邊為主的網(wǎng)絡。同時，為保證腦網(wǎng)絡的連通性，認為腦網(wǎng)絡中沒有孤立的節(jié)點[14]。因此，本文構建腦網(wǎng)絡時，所設閾值對應的鄰接矩陣的平均度(K)必須大于2ln90且網(wǎng)絡密度在10%～50%之間[14]。其中，網(wǎng)絡密度m為：

(2)

式中：N為網(wǎng)絡節(jié)點數(shù)；Ki為節(jié)點i的度數(shù)。根據(jù)以上規(guī)則，可得到描述腦網(wǎng)絡的鄰接矩陣。接著，再根據(jù)CLM相繼故障模型，對所構建的腦網(wǎng)絡進行模擬攻擊。

1.5　CLM相繼故障模型

CLM相繼故障模型由一般的無向無權圖G來表示傳輸網(wǎng)絡，各個節(jié)點i、j之間的邊權值eij∈[0,1]。eij表示該條邊上的傳遞效率。N×N的關聯(lián)矩陣{dij}代表了最短路徑值矩陣；N×N的關聯(lián)矩陣{eij}代表了節(jié)點之間的效率值矩陣。節(jié)點t時刻的負荷定義為：該時刻經(jīng)過該節(jié)點效率(距離)最優(yōu)(短)路徑的條數(shù)(該模型默認信息通過最短路徑傳遞)，其中負荷用Li(0)表示。效率最優(yōu)路徑定義為，對節(jié)點對(i，j)之間的所有路徑，計算該時刻效率值之和的最大值。同理，i、j兩點間最短路徑定義為，該時刻距離之和最小值。效率定義為：εij=1/dij，節(jié)點負荷容許值為：Ci=a×Li(0)(a≥1)，其中a代表容量參數(shù)。該模型用平均效率來表示整個系統(tǒng)的傳輸機制，即：

(3)

假設在t=0時刻，網(wǎng)絡G的傳輸效率值矩陣為{eij}，每個節(jié)點均有一個初始負荷。當對某個節(jié)點進行攻擊(在該模型下，攻擊某個節(jié)點亦為移除某個節(jié)點)時，必然會導致網(wǎng)絡G的效率最優(yōu)路徑改變，從而導致節(jié)點的負荷重新分配。可能引發(fā)節(jié)點過載，便會引起新一輪的負荷重新分配，最終可能導致相繼故障的發(fā)生。只有當每一次重新分配的負荷均小于對應節(jié)點的初始容量，該相繼故障才會停止。其中定義節(jié)點間的效率eij的演化公式為：

(4)

1.6　網(wǎng)絡特征分析

腦網(wǎng)絡經(jīng)CLM相繼故障模型攻擊后，其拓撲結構發(fā)生變化，全局效率與局部效率對網(wǎng)絡的拓撲結構具有重要描述意義。如下，重點描述全局效率與局部效率(整個網(wǎng)絡)以及某個節(jié)點的局部效率。全局效率定義為：

(5)

局部效率定義為：

(6)

式中：N表示網(wǎng)絡節(jié)點數(shù)，E(Gi)表示第i個節(jié)點的局部效率。E(Gi)定義為：第i個節(jié)點的所有鄰居節(jié)點間的效率均值。因此，整個網(wǎng)絡的局部效率為N個節(jié)點局部效率的均值。

2　結　果

2.1　構建腦網(wǎng)絡結果

本文閾值的設置從0.2開始，依次增加0.01，計算出不同閾值下對應的平均度。同時，必須滿足每個閾值對應的平均度(K)大于2ln90且不存在孤立節(jié)點。平均度與閾值關系如圖1所示，橫軸r表示閾值，縱軸K表示平均度值。隨著閾值增大，平均度逐漸減小，在閾值>0.31時，腦網(wǎng)絡開始出現(xiàn)孤立節(jié)點。同時，閾值為0.3時的平均度與網(wǎng)絡密度分別為14.69與16.5%，滿足平均度大于2ln90與網(wǎng)絡密度在10%～50%的條件。本文取閾值為0.3構造復雜腦網(wǎng)絡。由此，得到一個無權無相的鄰接矩陣，即表示構造的腦網(wǎng)絡模型。

圖1　閾值與平均度關系圖

2.2　負荷最大節(jié)點攻擊結果

負荷最大節(jié)點表示在t時刻最短路徑經(jīng)過該節(jié)點的次數(shù)最多。由于此節(jié)點的網(wǎng)絡拓撲較豐富，且具有高于全腦水平的局部效率，該節(jié)點對于整個網(wǎng)絡具有重要意義。本文構建的腦網(wǎng)絡中，負荷最大節(jié)點對應的腦區(qū)為左側眶內(nèi)額上回。

在不同容量系數(shù)下，對負荷最大的腦區(qū)進行攻擊，分析攻擊該節(jié)點之后的平均效率、整個網(wǎng)絡的局部效率、負荷最大節(jié)點的局部效率，得到圖2中所示的三條曲線。

圖2　不同容量系數(shù)下攻擊得到的全局效率、局部效率(整體、單個節(jié)點)變化趨勢

在圖2中，橫軸的a代表閾值，縱軸的E、表示效率值。其中，Eglob、Eloc、Eloc(nod)分別表示全局效率、整個網(wǎng)絡的局部效率、攻擊之后負荷最大節(jié)點的局部效率。三條曲線隨著容量系數(shù)a的增加幾乎均逐漸增加，且均逐漸逼近各自的初始效率(攻擊該節(jié)點后，網(wǎng)絡拓撲結構沒有發(fā)生改變時的效率)。同時，發(fā)現(xiàn)攻擊之后得到的全局效率、整個網(wǎng)絡的局部效率與負荷最大節(jié)點的局部效率在不同容量系數(shù)下依然具有較大值。在不同容量系數(shù)下，整個網(wǎng)絡的局部效率值高于全局效率值。且負荷最大節(jié)點的局部效率值在不同容量系數(shù)下高于全局效率值但小于整個網(wǎng)絡的局部效率值。

全局效率聚焦于整個網(wǎng)絡拓撲結構的變化，圖2中全局效率的變化趨勢顯示，容量系數(shù)與全局效率呈正相關。隨著容量系數(shù)遞增，腦網(wǎng)絡的全局效率越接近初始全局效率(0.503 1)。經(jīng)計算，當容量系數(shù)為1時，該腦網(wǎng)絡的全局效率較初始全局效率值降低52.6%。同時，發(fā)現(xiàn)容量系數(shù)為1.06時，該大腦全局效率較初始全局率值僅降低了8%，同比容量系數(shù)為1的大腦遭遇攻擊時的全局效率竟提高了93.7%。

3　討　論

攻擊結束之后，具有較高的全局效率與局部效率表明，整個腦網(wǎng)絡節(jié)點之間的信息交互在整個網(wǎng)絡與局部間具有較好的平衡。同時表明，腦網(wǎng)絡是一個具有較穩(wěn)定拓撲結構以及高效性能的網(wǎng)絡。該結論與Salvador等[15]的研究結果吻合。攻擊之后，整個網(wǎng)絡具有較高的局部效率表明正常人的腦網(wǎng)絡具有較好的容錯性[16]。結果中，全局效率在不同容量系數(shù)下均具有較高的效率值，表明腦網(wǎng)絡在全局網(wǎng)絡層面表現(xiàn)出較強的魯棒性。綜上所述，腦網(wǎng)絡對于抵抗相繼故障表現(xiàn)出較強的魯棒性且該網(wǎng)絡具有較穩(wěn)定的拓撲結構。對于整個腦網(wǎng)絡的局部效率高于全局效率的結果，表明局部節(jié)點間的路徑(效率)較短(較優(yōu))，在局部節(jié)點間的效率傳遞自然較高。但整個網(wǎng)絡間相對于局部節(jié)點間的路徑就變得較長，造成全局效率低于局部效率。進一步證明了腦網(wǎng)絡確實是一個稀疏且以大量短程邊為主的網(wǎng)絡[6]。同時，結果顯示攻擊結束后，新的負荷最大節(jié)點的局部效率低于整個網(wǎng)絡的局部效率。與初始網(wǎng)絡中負荷最大節(jié)點的局部效率高于整個網(wǎng)絡的局部效率截然相反。由此表明，攻擊該節(jié)點對于腦網(wǎng)絡的正常拓撲結構確實產(chǎn)生了一定的破壞，使腦網(wǎng)絡拓撲結構產(chǎn)生了較大的變化。提示，腦網(wǎng)絡雖具有較強魯棒性，但需注意保護腦網(wǎng)絡中負荷最大的節(jié)點，以維持較正常的拓撲結構。

此外，結果顯示，腦網(wǎng)絡中容量系數(shù)對于提升全局效率的效果顯著。表明，容量系數(shù)增強腦網(wǎng)絡的魯棒性。該結果與Kinney等[12]研究北美電力網(wǎng)絡時，發(fā)現(xiàn)容量系數(shù)對于平均效率的作用趨于一致。同時，Crucitti等[11]認為，容量系數(shù)相當于對一個復雜網(wǎng)絡進行“保護”，這種保護機制可以看作是相繼故障被系統(tǒng)本身“容納”和“吸收”了。從該相繼故障模型的角度，進一步分析容量系數(shù)的這種“容納”和“保護”機制。發(fā)現(xiàn)，每個節(jié)點均對應一個負荷，再乘以容量系數(shù)即是初始容量。因此，增加容量系數(shù)亦增加了該網(wǎng)絡的容負荷能力。而負荷是由經(jīng)過該節(jié)點的所有最短路徑次數(shù)決定。顯然，增大容量系數(shù)相當于縮短了該節(jié)點到達其他節(jié)點間的最短距離。這種距離的縮短可通過增加網(wǎng)絡冗余容量來實現(xiàn)。由此，增大容量系數(shù)，即增加網(wǎng)絡冗余容量，從而增強網(wǎng)絡的魯棒性。綜上可知，增加網(wǎng)絡冗余容量是一種抵抗相繼故障的有效手段。該結論與徐野等[17]發(fā)現(xiàn)小世界網(wǎng)絡的網(wǎng)絡冗余越大，其魯棒性越強的結論基本一致。因此，間接增加的網(wǎng)絡冗余容量才是“容納”和“吸收”機制的根源。不失一般性，增加腦區(qū)之間的冗余容量，復雜腦網(wǎng)絡的魯棒性必然增強。由腦區(qū)間功能連接強度來構造網(wǎng)絡可知，網(wǎng)絡冗余與功能連接強度成呈正相關。即增強某些腦區(qū)間的連接強度必然增強腦網(wǎng)絡的魯棒性。該研究結果可為腦網(wǎng)絡連接異常的疾病治療提供參考。

4　結　語

本文首次將CLM相繼故障模型應用到正常人腦網(wǎng)絡，該模型與腦網(wǎng)絡降低腦區(qū)間效率的方式較符合。利用Pearson相關系數(shù)構建腦網(wǎng)絡，分析它的功能連接水平并對其進行攻擊。本文發(fā)現(xiàn)，基于負荷最大攻擊策略下，腦網(wǎng)絡的魯棒性隨著容量系數(shù)的增大而增強。

對腦網(wǎng)絡的網(wǎng)絡拓撲屬性進行分析，發(fā)現(xiàn)腦網(wǎng)絡冗余越多，腦網(wǎng)絡魯棒性越強。提示增強腦區(qū)間的功能連接強度，可增強整個腦網(wǎng)絡的魯棒性。該結果為相關腦連接異常的疾病提供了治療參考。同時，發(fā)現(xiàn)整個腦網(wǎng)絡具有穩(wěn)定的拓撲結構以及較強魯棒性。并且提示注意保護腦網(wǎng)絡中負荷最大節(jié)點。

本文也存在諸多不足。比如，正常人樣本量較少。構建復雜腦網(wǎng)絡時，90個節(jié)點過于粗糙。未來的研究可從細分腦網(wǎng)絡節(jié)點著手研究。同時，本文構建的復雜腦網(wǎng)絡為無相無權圖。未來的研究可構建無相有權圖，使功能連接的劃分層次更精細。CLM相繼故障模型，在正常人腦網(wǎng)絡上的研究還需進一步推進。

[1] Watts D J,Strogatz S H.Collective dynamics of ‘small-world’ networks[J].Nature,1998,393(6684):440-442.

[2] Barabasi A L,Albert R.Emergence of scaling in random networks[J].Science,1999,286(5439):509-512.

[3] Stam C J.Functional connectivity patterns of human magnetoencephalographic recordings:a ‘small-world’ network?[J].Neuroscience Letters,2004,355(1):25-28.

[4] Salvador R,Suckling J,Schwarzbauer C,et al.Undirected graphs of frequency-dependent functional connectivity in whole brain networks[J].Philosophical Transactions of the Royal Society B Biological Sciences,2005,360(1457):937-946.

[5] Achard S,Salvador R,Whitcher B,et al.A resilient,low-frequency,small-world human brain functional network with highly connected association cortical hubs[J].Journal of Neuroscience,2006,26(1):63-72.

[6] Achard S,Bullmore E.Efficiency and cost of economical brain functional networks[J].PLoS Comput Biol,2007,3(2):17.

[7] Bassett D S,Bullmore E.Small-world brain networks[J].Neuroscientist,2006,12(6):512-523.

[8] 王艷群,李海芳,郭浩,等.靜息態(tài)腦功能網(wǎng)絡的社團結構研究[J].計算機應用,2012,32(7):2044-2048.

[9] 李越,郭浩,陳俊杰,等.抑郁癥功能腦網(wǎng)絡社團結構差異分析研究[J].計算機應用與軟件,2013,30(7):46-50.

[10] Pedersen M,Omidvarnia A H,Walz J M,et al.Increased segregation of brain networks in focal epilepsy:An fMRI graph theory finding[J].Neuroimage Clinical,2014,8(1):536-542.

[11] Crucitti P,Latora V,Marchiori M.Model for cascading failures in complex networks[J].Physical Review E Statistical Nonlinear & Soft Matter Physics,2004,69(4 Pt 2):045104.

[12] Kinney R,Crucitti P,Albert R,et al.Modeling cascading failur es in the North American power grid[J].The European Physical Journal B,2005,46(1):101-107.

[13] Song X W,Dong Z Y,Long X Y,et al.REST:A Toolkit for Resting-State Functional Magnetic Resonance Imaging Data Processing[J].Plos One,2011,6(9):e25031.

[14] 薛紹偉,唐一源,李健,等.一種基于fMRI數(shù)據(jù)的腦功能網(wǎng)絡構建方法[J].計算機應用研究,2010,27(11):4055-4057.

[15] Salvador R,Suckling J,Coleman M R,et al.Neurophysiological architecture of functional magnetic resonance images of human brain[J].Cerebral Cortex,2005,15(9):1332-1342.

[16] Sivan E,Parnas H,Dolev D.Fault tolerance in the cardiac ganglion of the lobster[J].Biological Cybernetics,1999,81(1):11-23.

[17] 徐野,王瑤.復雜網(wǎng)絡相繼故障的節(jié)點動態(tài)分析[J].沈陽理工大學學報,2015,34(1):17-21.