肖繼海，崔曉紅，陳俊杰

(太原理工大學信息與計算機學院，山西 晉中 030600)

1 引言

近年來，圖挖掘已經(jīng)成為一個熱門的研究領域并應用到社交網(wǎng)絡、計算生物學及計算機網(wǎng)絡的分析中。此外，許多新的數(shù)據(jù)也可以表示為圖，如功能性磁共振大腦影像數(shù)據(jù)fMRI(functional Magnetic Resonance Imaging)。本文利用受試者的fMRI數(shù)據(jù)，構建腦網(wǎng)絡，其中每個節(jié)點代表一個腦區(qū)，每條邊代表2個腦區(qū)之間的功能連通性[1]。這些腦網(wǎng)絡為探索人腦的內(nèi)部結構和活動提供了一個新的視角,同時也可為腦疾病臨床診斷提供有價值的輔助信息。

目前的研究大都集中于腦疾病的分類，例如，Huang等人[2]提出了一種基于樹引導的組稀疏學習多頻帶腦網(wǎng)絡分類模型，首先在多頻帶中獲取低頻振幅比率數(shù)據(jù)，使用基于樹的稀疏學習方法進行特征選擇，最后結合多頻段特征進行分類。Liu等人[3]從T1加權MRI(Magnetic Resonance Imaging)圖像和DTI(Diffusion Tensor Imaging)圖像中提取8個特征集，并采用兩步特征選擇方法得到最具鑒別性的特征，最后采用基于多特征集的多核SVM學習方法對精神分裂癥SZ(SchiZophrenia)患者和正常人NC(Normal Control)進行分類。楊勇哲等人[4]提取結構和靜息態(tài)功能磁共振圖像的90個感興趣區(qū)域的灰質體積、局部一致性、低頻振蕩振幅和度中心性作為特征，使用基于遞歸特征消除的支持向量機對腦疾病患者和正常人進行分類。

以上研究旨在選擇特征來訓練分類器，以實現(xiàn)對數(shù)據(jù)的自動分類。當訓練數(shù)據(jù)樣本量足夠大并且選擇了最優(yōu)特征時，分類性能通常很好。但是，從大量的數(shù)據(jù)中提取和選擇最優(yōu)的特征很費時。此外，在許多實際情況下，只能獲取到較少的標記數(shù)據(jù)樣本。因此，在有限的標記數(shù)據(jù)樣本上找到有顯著差異的特征是不可靠的。

針對以上問題，本文提出了一種節(jié)點屬性和拓撲信息相結合的腦網(wǎng)絡聚類模型。與分類相反，聚類的不同之處在于不需要數(shù)據(jù)的類型標簽[5]。腦網(wǎng)絡聚類的目的是根據(jù)腦網(wǎng)絡之間的相似程度自動聚類，類內(nèi)腦網(wǎng)絡之間具有較高的相似性[6]。

該方法主要的挑戰(zhàn)是如何準確計算腦網(wǎng)絡之間的相似度。目前，文獻[7]調查了48名SZ患者和24名健康對照者的全局和節(jié)點網(wǎng)絡特性，發(fā)現(xiàn)精神分裂癥患者的大腦網(wǎng)絡特性與健康人的相比有顯著差異，比如SZ患者的平均介數(shù)中心性較高。另外，文獻[8]發(fā)現(xiàn)SZ患者功能性磁共振網(wǎng)絡的拓撲和空間結構被顯著破壞。以上文獻表明，SZ患者的腦網(wǎng)絡在節(jié)點屬性和拓撲結構上都存在顯著差異。因此，為了更全面地捕獲腦網(wǎng)絡的差異，本文從節(jié)點屬性和拓撲結構兩方面測量腦網(wǎng)絡的相似度。首先，選擇目前研究公認的SZ患者存在異常功能連接的默認模式網(wǎng)絡DMN(Default Mode Network)[9,10]作為感興趣區(qū)域構建腦網(wǎng)絡，然后分別使用余弦相似度和子網(wǎng)絡核來度量腦網(wǎng)絡的屬性相似度和結構相似度，接著將結構相似度和屬性相似度集成為一個相似度矩陣，最后利用譜聚類實現(xiàn)腦網(wǎng)絡聚類。

2 腦網(wǎng)絡聚類模型

腦網(wǎng)絡聚類模型如圖1所示，包含fMRI數(shù)據(jù)的預處理、標記DMN、構建相應腦網(wǎng)絡、計算腦網(wǎng)絡的相似度以及聚類。也就是，逐一對每個腦網(wǎng)絡的fMRI數(shù)據(jù)先實施預處理，找到DMN對應的腦區(qū)并構造腦網(wǎng)絡，隨后計算其屬性相似度和結構相似度，最后整合2種不同類型的相似度，并使用譜聚類方法對腦網(wǎng)絡進行聚類。

Figure 1 Brain network clustering model圖1 腦網(wǎng)絡聚類模型

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

使用DPARSF(Data Processing Assistant for Resting-State fMRI)工具包對fMRI圖像數(shù)據(jù)進行預處理，包括刪除前10個時間點數(shù)據(jù)、時間片的校正和頭動的校正。本文使用DARTEL(Diffeomorphic Anatomical Registration Through Exponentiated Lie algebra)配準方法將數(shù)據(jù)配準到標準空間，期間將數(shù)據(jù)重新采樣為3×3×3 mm3體素，再進行平滑處理以去除協(xié)變量的影響，最后，進行濾波去噪，頻率范圍為[0.01,0.08] Hz，以減少高頻和極低頻生理噪聲的影響，這樣能較好地反映神經(jīng)元自發(fā)活動。

2.2 標記DMN并構建功能腦網(wǎng)絡

研究已證實，在靜息態(tài)功能腦網(wǎng)絡中默認模式網(wǎng)絡(DMN)能保持相對穩(wěn)定的狀態(tài)，通過它可以研究功能腦網(wǎng)絡連接的異常[10 - 12]。本文使用自動化解剖標記AAL(Automated Anatomical Labeling)圖譜[13]，從中找到DMN所包含的32個腦區(qū)(左右腦)，腦區(qū)名稱如表1所示。

Table 1 Brain regions included in the DMN

腦網(wǎng)絡的節(jié)點由表1中的32個腦區(qū)定義，節(jié)點的時間序列取腦區(qū)內(nèi)所有體素時間序列的平均值，節(jié)點之間的連接權值由節(jié)點時間序列間的皮爾遜相關系數(shù)定義，加權全連接網(wǎng)絡由此形成。為了避免傳統(tǒng)二值網(wǎng)絡構建過程中產(chǎn)生的偏差，使用Kruskal算法[14]將加權全連接網(wǎng)絡轉換成無偏差的腦網(wǎng)絡，最終形成最小生成樹無偏DMN腦網(wǎng)絡。

2.3 測量腦網(wǎng)絡相似度

腦網(wǎng)絡不僅具有全局拓撲結構特征而且還具有局部屬性特征。為了綜合評價腦網(wǎng)絡的相似度，本文將上述2方面特征結合起來測量腦網(wǎng)絡相似度，流程圖如圖2所示。圖2中，G和H為腦網(wǎng)絡，S(G,H)為腦網(wǎng)絡G和H的相似度，Satt(G,H)為G和H的屬性相似度，Sstr(G,H)為G和H的結構相似度，δ為權重。

Figure 2 Flow chart of similarity measurement in brain network圖2 腦網(wǎng)絡相似度測量流程圖

(1)節(jié)點屬性相似度的測量。

介數(shù)屬性是復雜網(wǎng)絡的一個重要特征。節(jié)點介數(shù)越大，表示該節(jié)點是中樞節(jié)點，對整個網(wǎng)絡的信息傳輸起到很大的作用。研究表明[15]，SZ患者的大腦中樞受到破壞，在臨床應用中，介數(shù)常用于健康人和SZ患者腦網(wǎng)絡的比較。因此，本文使用節(jié)點的介數(shù)來表征腦網(wǎng)絡的局部屬性，腦網(wǎng)絡間屬性相似度由腦網(wǎng)絡中所有節(jié)點介數(shù)之間的相似度來評價。通過節(jié)點i的所有最短路徑的數(shù)量定義為節(jié)點i的介數(shù)bci，計算如式(1)所示：

(1)

屬性相似度Satt(G,H)使用余弦相似度計算，如式(2)所示：

(2)

其中,bcm(H)表示腦網(wǎng)絡H中第m個節(jié)點的介數(shù)，bcm(G)表示腦網(wǎng)絡G中第m個節(jié)點的介數(shù)，n表示腦網(wǎng)絡中節(jié)點總數(shù)。

(2)拓撲結構相似度測量。

一組對象之間的相似性可以用核函數(shù)來度量。當用核函數(shù)處理圖數(shù)據(jù)時，該核稱為圖核，圖核是把圖數(shù)據(jù)從原始的圖空間映射到特征空間，使得2幅圖像的相似性就是它們在特征空間中的點積。因此，2幅圖之間的相似性可通過圖核來度量[16]。

大腦網(wǎng)絡與圖的主要區(qū)別是大腦網(wǎng)絡中每個節(jié)點的唯一性，也就是說，大腦網(wǎng)絡中的每個節(jié)點都表示特定的大腦區(qū)域，是大腦網(wǎng)絡的一個特定特征。因此，本文采用子網(wǎng)絡核(Sub-network Kernel)[17]測量腦網(wǎng)絡間的拓撲結構相似度。與傳統(tǒng)意義上的圖核相比，子網(wǎng)絡核的一個重要特性就是能捕獲網(wǎng)絡節(jié)點上的多層次拓撲性質，而且還兼顧了網(wǎng)絡中各個節(jié)點的唯一性。子網(wǎng)絡核的構建過程簡要描述如下：

①在每個網(wǎng)絡節(jié)點上構建一組子網(wǎng)絡，以便反映腦網(wǎng)絡在多層次上的連通性。

G=(V,E)和與H=(V,E′)表示一對腦網(wǎng)絡，其中V代表腦網(wǎng)絡的節(jié)點集，E和E′分別表示G與H的邊集。各個腦網(wǎng)絡擁有相同的大腦區(qū)域，即它們有相同的節(jié)點。為了表達腦網(wǎng)絡的多級拓撲特性，需要在網(wǎng)絡G與網(wǎng)絡H的每個節(jié)點vi上分別定義2組子網(wǎng)絡：

(3)

(4)

由式(3)和式(4)可知，對于由n個節(jié)點組成的腦網(wǎng)絡，可以得到n組子網(wǎng)絡：

(5)

(6)

②網(wǎng)絡G與H中同一節(jié)點的相似度就是2個網(wǎng)絡中同一節(jié)點的所有子網(wǎng)絡的相似度，如式(7)所示：

(7)

(8)

(9)

(10)

③因此，腦網(wǎng)絡G與H的圖核就是腦網(wǎng)絡中所有節(jié)點對應子網(wǎng)絡對的相似度，n表示腦網(wǎng)絡中的節(jié)點數(shù)目，G和H的圖核如式(11)所示：

(11)

④最后得到2個腦網(wǎng)絡G與H的圖核即為2個腦網(wǎng)絡G與H之間的拓撲結構相似度，如式(12)所示：

Sstr(G,H)=k(G,H)

(12)

(3)腦網(wǎng)絡的相似度。

為了全面評價腦網(wǎng)絡之間的相似度，本文將腦網(wǎng)絡之間屬性方面的相似度與結構方面的相似度進行優(yōu)化組合。

首先，因為屬性相似度與結構相似度的類型不同，故在組合二者之前必須歸一化，如式(13)所示：

(13)

其中，Snorm(G,H)是腦網(wǎng)絡G與H歸一化之后的相似度，S(G,H)是腦網(wǎng)絡G與H的相似度。

然后，通過權重δ來調整2種相似度的貢獻程度，使組合最優(yōu)。

腦網(wǎng)絡G與H的相似度定義如公式(14)所示：

SG,H=δSatt_norm(G,H)+

(1-δ)Sstr_norm(G,H)

(14)

其中，Satt_norm(G,H)是功能腦網(wǎng)絡G與H歸一化后的屬性相似度，Sstr_norm(G,H)是功能腦網(wǎng)絡G與H歸一化后的結構相似度，δ是權重參數(shù)。

第三步，腦網(wǎng)絡之間的相似矩陣S定義如下：

(15)

其中,SG,H是腦網(wǎng)絡G與H之間的相似度，N是腦網(wǎng)絡數(shù)量。

2.4 腦網(wǎng)絡聚類模型描述

前文已經(jīng)構造出腦網(wǎng)絡的相似矩陣S，這樣可以方便地將腦網(wǎng)絡的聚類問題轉化成譜聚類[19]問題，將相似度高的腦網(wǎng)絡自動歸為一類。聚類模型的基本過程描述如下：

輸入：N個腦網(wǎng)絡。

輸出：向量Y，若向量Y的第i個元素為cj，則表示第i個被試對應的腦網(wǎng)絡歸屬到類cj中。

步驟1將腦網(wǎng)絡的相似矩陣S置為N×N零元矩陣，N表示腦網(wǎng)絡數(shù)目；

步驟2對N個腦網(wǎng)絡依次計算每個節(jié)點的介數(shù)；

步驟3由式(2)得到N個腦網(wǎng)絡之間的屬性相似度；

步驟4對N個腦網(wǎng)絡執(zhí)行子網(wǎng)絡核算法，得到N個腦網(wǎng)絡之間的結構相似度；

步驟5將屬性相似度和結構相似度歸一化，由式(14)與式(15)得到相似矩陣S；

步驟6定義對角陣D，其元素D(i,i)為相似矩陣S的第i行元素之和；

步驟7定義矩陣L=D-1/2SD-1/2;

步驟8求矩陣L的m個最大特征向量x1,x2,…,xm，按列存放m個特征向量產(chǎn)生矩陣X=[x1,x2,…,xm]∈Rn×m；

步驟9將矩陣X的每一行標準化，得到矩陣Y；

步驟10將矩陣Y的每一行作為Rm中的一個點，執(zhí)行K-means算法將這些點聚成m個類。

3 實驗數(shù)據(jù)與結果

3.1 實驗數(shù)據(jù)

本文實驗的數(shù)據(jù)來自于openfMRI數(shù)據(jù)庫(https://www.openfmri.org/)，經(jīng)篩選，滿足數(shù)據(jù)收集條件的有99名，其中包含50名SZ患者與49名正常對照組。掃描參數(shù)：TR(Repetition Time)重復時間為2 000 ms；TE(Echo Time)回波時間為30 ms；翻轉角度為90°；切片厚度為4 mm；切片數(shù)量為34；時間點數(shù)為152個。被試的基本信息統(tǒng)計情況如表2所示。

Table 2 Statistical table of basic information of subjects

3.2 實驗參數(shù)設置

實驗中，因為僅有2個可能的標簽(SZ，NC)，故將類別數(shù)m設置為2。在子網(wǎng)絡核算法中，參數(shù)t與d分別設置為3與3。權重δ通過網(wǎng)格搜索法來優(yōu)化，δ初始值是0.1，步長是0.1，最大值是0.9。

3.3 無偏腦網(wǎng)絡

為了構建合理的腦網(wǎng)絡，本文腦網(wǎng)絡的節(jié)點由DMN包含的32個腦區(qū)定義，節(jié)點之間的連接權值由節(jié)點時間序列間的皮爾遜相關系數(shù)定義。使用Kruskal算法得到SZ與NC所有被試的無偏腦網(wǎng)絡。SZ(No.01)與NC(No.01)的無偏腦網(wǎng)絡如圖3所示。從圖3中看出，SZ患者和正常人的大腦在額葉、顳葉和頂下小葉之間的功能連接存在明顯的差異。

Figure 3 Unbiased brain network of SZ and NC圖3 SZ與NC的無偏腦網(wǎng)絡

3.4 相似矩陣

譜聚類的關鍵在于相似性度量。腦網(wǎng)絡構建后，腦網(wǎng)絡之間的相似度通過整合腦網(wǎng)絡在屬性與結構2方面的相似度來綜合評價。本文首先使用余弦相似度來計算腦網(wǎng)絡間的介數(shù)屬性相似度，腦網(wǎng)絡的屬性相似度矩陣如圖4所示；其次，使用子網(wǎng)絡核來度量腦網(wǎng)絡間的結構相似性，其相似度矩陣如圖5所示；最后，選擇最優(yōu)權重參數(shù)δ，將屬性相似度與結構相似度進行優(yōu)化組合，得到最終相似度矩陣，如圖6所示。

Figure 4 Attribute similarity matrix圖4 屬性相似度矩陣

3.5 聚類性能

將屬性相似度矩陣與結構相似度矩陣按一定比例優(yōu)化組合，形成最終相似度矩陣。然后通過譜聚類[20]完成上述腦網(wǎng)絡的聚類。通過以下指標來評估聚類方法的性能： Rand指數(shù)、召回率、精確率和F-measure[21]。為準確評價上述指標，實驗中要求重復運行程序50次，將其平均值作為最終的性能評價結果。

Figure 5 Structural similarity matrix圖5 結構相似度矩陣

Figure 6 Final similarity matrix圖6 最終相似度矩陣

(1)Rand指數(shù)表示正確聚類的百分比，Rand指數(shù)越大，意味著聚類結果與真實情況越吻合。

Rand指數(shù)由式(16)計算：

(16)

其中,TP是真正類的個數(shù)，TN是真負類的個數(shù)，F(xiàn)P是假正類的個數(shù)，F(xiàn)N是假負類的個數(shù)。

(2)召回率(Recall)表示在被識別為正類別的樣本中，真實情況是正類別的比例，即查全率，由式(17)計算：

(17)

(3)精確率(Precision)表示在所有正類別樣本中，被正確識別為正類別的比例，即查準率,由式(18)計算：

(18)

(4)F-measure綜合了精確率和召回率的結果，是精確率與召回率加權調和平均值。F-measure由式(19)定義：

(19)

其中,P表示精確率，R表示召回率。在本文中，使用F1(β=1)評價聚類的性能。當F1的值較高時，則說明聚類模型的效果比較理想。

為評估本文所提方法SA-cluster(Structure Attribute-cluster)的聚類性能，在相同數(shù)據(jù)集上設計實驗,將本文的方法與2種方法進行比較：

(1)基于介數(shù)屬性的譜聚類算法(A-cluster)：構造基于介數(shù)屬性的相似矩陣，將相似矩陣輸入譜聚類模型，對腦網(wǎng)絡進行聚類。

(2)基于結構的譜聚類方法(S-cluster)：構造基于拓撲結構的相似矩陣，將相似矩陣輸入譜聚類模型，對腦網(wǎng)絡進行聚類。

相同數(shù)據(jù)集上3種方法的聚類性能如表3所示。結果表明，本文方法(SA-cluster)的聚類性能最好，其Rand指數(shù)是0.91，精確率是0.86，召回率是0.98，F(xiàn)1是0.92。

Table 3 Clustering performance of three methods

4 實驗結果討論

4.1 聚類性能的比較

如表3所示，本文方法在Rand指數(shù)、精確率、召回率和F1方面表現(xiàn)最出色。3種方法中，A-cluster方法使用腦網(wǎng)絡的介數(shù)屬性計算每對腦網(wǎng)絡之間的相似性，并用于聚類。A-cluster方法的Rand指數(shù)與召回率比較低，這說明該方法不能較準確地識別出SZ患者，存在漏診。結果表明，僅僅從局部屬性角度不能準確反映腦網(wǎng)絡之間的相似性。

S-cluster通過在腦網(wǎng)絡上執(zhí)行子網(wǎng)絡核算法來得到相似矩陣。從表3可以看出，與A-cluster方法相比，該方法的聚類結果有所改善，查準率和查全率提高了，漏診率下降，能較準確地識別出患者。結果說明，當網(wǎng)絡的局部屬性發(fā)生變化時，網(wǎng)絡的全局連接結構也受到影響。因此，全局結構的相似性描述方法對聚類有很大貢獻。

與前2種方法相比，本文方法(SA-cluster)的Rand指數(shù)是0.91，精確率是0.86，召回率是0.98，F(xiàn)1是0.92，準確率和查全率最高，漏診率最低，更能準確地識別出患者，聚類的效果比較理想。這主要是因為節(jié)點屬性和拓撲結構從不同的角度描述腦網(wǎng)絡的相似度，將兩者有效地結合可以更好地表達腦網(wǎng)絡之間的相似度。這一結論正好與文獻[7,8]的研究結果一致，進一步證實了SZ患者的大腦網(wǎng)絡在節(jié)點屬性和拓撲結構上與正常人的腦網(wǎng)絡存在差異。

基于上述分析，本文方法巧妙地整合了屬性相似度與結構相似度，其中以介數(shù)屬性相似度來刻畫腦網(wǎng)絡的局部特征，以結構相似度來刻畫腦網(wǎng)絡的全局拓撲結構特征。特別地，本文使用子網(wǎng)絡核算法來度量腦網(wǎng)絡的拓撲結構相似性，該算法不僅考慮到腦網(wǎng)絡中每個腦區(qū)節(jié)點的唯一性，而且還能夠捕獲到腦區(qū)節(jié)點的多級拓撲特征，這一點對于腦網(wǎng)絡的相似性度量非常重要。可以看出，屬性特征與內(nèi)部節(jié)點拓撲結構對圖聚類具有重要意義。因此，基于屬性與結構相結合的相似性度量方法可以更準確地刻畫腦網(wǎng)絡的相似性，從而有效提高聚類性能。

4.2 參數(shù)δ對聚類性能的影響

權重參數(shù)δ用于調整屬性相似度與結構相似度的貢獻程度。為分析權重參數(shù)δ對本文方法聚類性能的影響，本文為參數(shù)δ設置了不同的值，δ初始值是0.1，遞增步長是0.1，終值是0.9。參數(shù)δ取不同值的聚類結果如圖7所示。

Figure 7 Influence of different δ values on clustering performance圖7 參數(shù)δ不同取值對聚類性能的影響

圖7表明，無論δ取何值，Rand指數(shù)、精確率、召回率和F1值都在0.8以上，說明從屬性和結構2方面能較準確地衡量腦網(wǎng)絡相似度，腦網(wǎng)絡的拓撲結構和屬性特征對腦網(wǎng)絡的聚類起著重要作用。當δ取值為0.2時，聚類性能最好，其Rand指數(shù)為0.91,精確率為0.86,召回率為0.98,F1值為0.92，這說明在相似度的計算中，結構相似度貢獻較大。

5 結束語

本文提出一種節(jié)點屬性和拓撲信息相結合的腦網(wǎng)絡聚類模型。首先，使用余弦相似度來計算腦網(wǎng)絡之間的屬性相似度，然后，使用子網(wǎng)絡核來度量腦網(wǎng)絡之間的結構相似度，最后，將屬性相似度與結構相似度相結合，得到相似度矩陣，并用譜聚類算法聚類。實驗結果表明，本文提出的聚類模型在Rand指數(shù)、召回率、精確率和F1值方面都表現(xiàn)較優(yōu)。