熵模型廣泛應用在癲癇腦電分析中，但其在靜息態功能磁共振成像（rfMRI）中尚存在主觀選擇計算參數的問題。為此，本文提出多尺度熵模型優化算法，聯合有監督機器學習檢驗優化效果。以致癇側定位為例，將20位海馬硬化標記患者分為左、右側2組，利用敏感性分析指標優化熵模型參數，以組間優化熵值有顯著差異的腦區作為對致癇側敏感的標記，其熵值為特征向量輸入支持向量機分類并驗證，獲得平均準確率達95%的定側結果，高于目前水平。研究結果顯示，熵模型參數優化算法可較為準確地提取對致癇側敏感的功能影像學標記，達到客觀選擇癲癇rfMRI熵模型參數的研究目的，為熵應用于先進技術檢測提供了依據。
引用本文： 謝沖, 葛曼玲, 付曉璇, 陳盛華, 張夫一, 郭志彤, 張志強. 靜息態功能磁共振成像多尺度熵優化及在顳葉癲癇致癇側定位中的應用. 生物醫學工程學雜志, 2021, 38(6): 1163-1172. doi: 10.7507/1001-5515.202011048 復制
引言
癲癇是由大腦神經元異常放電引起的中樞神經系統功能異常的慢性疾病，嚴重影響著患者生活和生命質量，且病變可在大規模腦網絡中產生分布效應，導致與原發病灶相連的遠端腦區發生代謝功能障礙或者結構改變[1]。手術治療是難治性癲癇的首選方法，70%以上的局灶性癲癇患者手術治療后病情得到緩解，顳葉癲癇（temporal lobe epilepsy，TLE）患者表現出最積極的治療效果[2-3]。因此，癲癇致癇側、致癇灶和癲癇傳播路徑等定位是術前的重要準備，定側是基礎[4]。
臨床上，對于癲癇疾病診斷以及致癇灶定位常采用腦電圖和常規結構磁共振成像（magnetic resonance imaging，MRI）技術。腦電圖是臨床上定位致癇灶和致癇側的金標準，其通過監測始發異常神經放電實現定位。然而，顱內電極或立體定位電極的插入會不可避免地造成腦組織損傷，并且腦電存在空間分辨率低、靈敏度不足、觀察范圍有限和費用昂貴等缺點[5]。此外，常規MRI作為TLE患者術前定位最常用的檢查方法，通過對海馬形態的觀察和海馬體體積測量能夠對部分TLE患者作出定位診斷[6]。海馬陽性是致癇側的典型標識，但只有海馬區神經元丟失達50%機體形態才會表現陽性的異常標記[7]，所以，常規MRI掃描對致癇灶和致癇側定位診斷的檢出率（33%）較低[8]。
近年來，靜息態功能磁共振成像（resting-state functional MRI，rfMRI）以其無損、快捷和廣泛的患者適應性等優勢，成為檢測癲癇引起異常腦功能的先進技術。研究表明，不同致癇側的TLE患者rfMRI存在較大差異[9-10]。Barron等[11]發現不同致癇側的TLE患者丘腦和顳葉之間功能連接強度存在差異；Yang等[12]分析TLE患者的靜息態功能連接，發現該方法選取的特征量能夠較為有效地區分TLE患者的致癇側。rfMRI不僅可以從腦功能水平揭示疾病腦區病理變化，還能增加手術的安全性和準確性，對神經疾病術后及療效監測有較大的臨床應用價值[13-15]。
以往，熵模型廣泛應用于癲癇腦電信號分析。為從神經信號中提取更多有效信息，Costa等[16]提出多尺度熵（multi-scale entropy，MSE）模型將熵和尺度聯合共同表征生理系統在不同狀態下所表現的復雜特性。MSE模型被廣泛用于分析具有多尺度和非線性特點的腦電[17]、fMRI[18-19]、功能性近紅外光譜（functional near-infrared spectroscopy，fNIRS）[20]和心腦信號[21-23]。Yang等[24]通過MSE模型發現年輕人和老年人的血氧水平依賴（blood oxygen level dependent，BOLD）信號復雜度存在差異，Li等[20]通過MSE模型發現輕度認知障礙患者、阿爾茨海默病患者和健康人的fNIRS信號復雜度也存在差異。然而，MSE模型的參數選擇對信號本身存在較大依賴性，以往多依賴經驗選取參數，具有一定的主觀性。對rfMRI信號來說，常采用m = 1、r = 0.5[19]或m = 1、r = 0.35[24]等。
TLE患者rfMRI信號的非線性顯著高于健康人[25-26]，因此，應用MSE模型研究TLE患者的rfMRI復雜度，或可為患者的檢測或治療提供新思路。本文以rfMRI致癇側定位為例，試圖用敏感性分析的接受者操作特征曲線（receiver operating characteristic，ROC）及其曲線下面積AUC值（area under ROC curve，AUC）優化MSE模型的計算參數。為檢驗優化效果，通過優化后的熵模型提取對致癇側敏感的功能像標記，并將其熵值作為特征向量輸入至可適應小樣本的機器學習模型——支持向量機（support vector machine，SVM）定側，以其精度評估優化有效性，試圖解決MSE參數選擇依賴經驗和數據特征、缺少統一客觀依據的計算瓶頸，為熵值應用于TLE患者rfMRI提供客觀依據。最后，由于SVM是有監督學習模型，需要在正確標注的樣本上運行，所以，本文借助結構像海馬陽性提示的致癇側信息進行有監督的機器學習。整體分析框架如圖1所示。
圖1
分析流程圖
Figure1.
Study flowchart
1 數據與方法
1.1 數據采集
本文中20名TLE患者rfMRI數據來自南京大學醫學院附屬金陵醫院醫學影像科，其中男性9例，女性11例，年齡19～33（25.7±3.7）歲，接受常規的術前評估以定位癲癇放電區域。其中結構像海馬陽性左側和右側的患者各10例（男性，左側4例，右側5例；女性，左側6例，右側5例），命名為左側組和右側組。
數據采集使用Siemens 3.0T高場強超導MR掃描儀，經單次激發GE-EPI序列采集BOLD-fMRI數據，TR = 2 000 ms，TE = 30 ms，翻轉角度為90°，矩陣64×64，掃描視野240 mm×240 mm，層厚4.00 mm，層間距0.40 mm[27]。本研究獲得南京軍區總醫院醫學倫理委員會批準，患者知情并且簽署知情同意書。
1.2 數據預處理
使用FMRIB軟件庫中的工具進行rfMRI數據預處理[28-29]，處理流程和質量檢測（頭動控制量 < 0.2 mm）在美國哈佛醫學院Martinos醫學影像中心完成[30-31]。處理平臺CPU為Intel Xeon Sliver 4112 × 16，64核，操作系統Centos7.6，每個被試者的數據預處理時間約為15 h。
預處理步驟：① 將每個被試采集的前4個時間點數據移除，使信號更加穩定；② 使用SPM2校正切片采集的每個腦體積時間偏移（Wellcome Department of Cognitive Neurology，London，英國）；③ 利用FSL軟件（http://fsl.fmrib.ox.ac.uk/fsl）對頭部運動進行剛體平移和旋轉[32]；④ 全腦平均信號的歸一化，將信號配準到蒙特利爾神經研究所（Montreal Neurological Institute，MNI）的標準空間，重采樣到2 mm各向同性體素大小；⑤ 帶通時間濾波（0.01～0.08 Hz）。
預處理后，使用win10操作系統MatlabR2018b版本的DPABI將數據投射到自動解剖標記（the anatomical automatic labeling atlas，AAL1）上，每個被試者的數據處理時間約為20 min，將每個被試全腦分成116個腦區，每個腦區有246個時間點，時間長度為490 s。
1.3 多尺度熵
1.3.1 多尺度熵模型
多尺度熵模型通過評估量化生物醫學信號的“復雜度”，能夠有效分析信號的動態變化，因此在神經信號分析中具有應用價值，特別是癲癇這種動態疾病。
多尺度熵的計算步驟：
1）將采集長度為N的時間序列
進行粗粒化，重新構建一個時間序列
如式（1）所示。粗粒化（非重疊式）是根據選取的時間尺度大小，對連續的時間點取平均值，重新構建時間序列。
 |
式中
為時間尺度。
2）定義多尺度熵模型參數m，其中m為整數，代表比較向量的長度。按順序構建一組m維矢量，即
，
。
3）定義多尺度熵模型參數公差閾值（即相似系數r），其中r為實數，代表匹配閾值。定義Ui與Uj之間的距離 
為兩個矢量中元素之間距離的最大值，即 
，其中i ≠ j，統計當 
時 
小于r時的數目b，并將b與距離總數 
的比值，記做 
，如式（2）所示。
 |
4）將維數m加1，在重新設定維度m的情況下重復步驟2、3，將會得到
。
5）理論上，多尺度熵計算如式（3）所示：
 |
當N有限時，式（3）的計算模型如式（4）所示：
 |
綜上所述，多尺度熵中的3個參數分別為維度m、相似度r和時間尺度τ。若計算所得熵值越大，表明復雜度越大，反之亦然。
1.3.2 熵模型參數優化準側
多尺度熵是在不同時間尺度值下計算時間序列的樣本熵，所以，多尺度熵模型中維度m和相似度r的選取與樣本熵參數選取相似。因為多尺度熵的3個參數之間相互影響，則在單個參數優化過程中，其余兩個參數也須參與優化過程，即3個參數同時參與一個參數的優化過程。本文中，首先利用組間熵值差異較大（P ≤ 0.05）的腦區數量確定維度m和另外兩個參數即將優化的區間，然后，結合ROC曲線和AUC值依次對相似度r和時間尺度 τ 進行優化。多尺度熵模型的參數優化流程見圖2。
圖2
模型參數優化流程圖
Figure2.
Flowchart of parameters optimization in a multi-scale entropy model
原則上，熵評估的準確度會隨著維度m和m + 1的矢量（時間序列）匹配數量的增加而提高，選擇較小的維度m和較大的相似度r也能夠增加矢量匹配數量提高準確度。基于時間長度10m～20m（m參數維度）的BOLD信號可較為準確地計算熵值[33]。根據本文中數據的時間點長度可以選取m = 1或2。
根據經驗，相似度r范圍可選在0.05～0.6。但對文中數據分析發現：在r = 0.05～0.3時，腦區熵值會出現無效值，這是由于選取相似度r較小，計算過程中不存在與之相匹配的矢量。為保證分析精確性，相似度r選在0.3～0.6（步長0.02）范圍內，時間尺度也不宜選擇過大，本文選取范圍是 
= 1～5。
本文用敏感性分析的ROC曲線（反映分類效果敏感性和特異性連續變量的綜合指標）和AUC值（ROC曲線下方的面積大小）來評估MSE模型參數的優化，使參數有利于分析癲癇rfMRI信號特征，使得參數選擇具有一定的統計特性和客觀性。
1.3.3 對發作側敏感的腦區
利用MSE模型計算全腦（AAL1）中90個大腦區域的rfMRI信號得到腦區熵值，并通過t檢驗對比兩組患者腦區熵值的差異性，獲得腦區P值。在P ≤ 0.05范圍內，選取腦區將其作為對致癇側敏感的功能影像標記，并利用ROC曲線和AUC值（IBM SPSS Statistics21，美國）對標記進行驗證。最后用BrainNet Viewer[34]（http://www.nitrc.org/projects/bnv/）將標記進行可視化。
1.4 機器學習
1.4.1 特征向量
在優化參數下，將對致癇側敏感的腦區熵值作為特征向量，輸入到支持向量機中。
1.4.2 支持向量機
由Vapnik提出的SVM是一種有監督的機器學習方法，以訓練誤差作為優化問題的約束條件，并以置信范圍最小化為優化目標，是一種基于結構風險最小化準則的學習方法，不易受數據維數和樣本數量的限制，這為其在小樣本中的應用提供了依據[35-36]。
文中將20名患者的腦區多尺度熵值作為特征向量，按照3∶2比例分成12名訓練集和8名測試集，作為機器學習的基本輸入信息進行數據歸一化。選擇徑向基函數（radial basis function，RBF）做為核函數，選取比例參數g = 
代替核參數
，形成一組參數對（C，g）。其中，參數C和g的取值范圍設定為［? 10，10］（步長0.2），并用網格搜索方法（交叉驗證法）尋找誤差最小的懲罰系數C和核參數
為最優參數，創建訓練模型，對測試集數據進行計算。
1.4.3 留一交叉驗證
留一交叉驗證（leave one out cross-validation，LOOCV）常被用來檢測機器學習模型的準確性，其樣本利用率較高，能夠從有限的樣本中盡可能地獲取更多有效信息，所以，適用于小樣本分析。如果有N個樣本，其中N ? 1個為訓練樣本，剩余1個為測試樣本，用N次結果的平均分類精度來評估模型的準確性（文中，N = 20）。
2 結果
2.1 MSE模型的優化參數
調節3個參數（維度m、相似度r和時間尺度τ）計算腦區熵值，通過t檢驗分析左、右致癇側TLE患者間的腦區熵值，從而獲得具有統計學意義（P ≤ 0.05）的腦區和優化參數范圍，并利用ROC曲線和AUC值對其驗證。
2.1.1 維度m優化值
以時間尺度τ = 1～5（步長1）生成圖，每張圖以相似度r = 0.3～0.6（步長為0.02）為變量，對比維度m = 1和m = 2時組間熵值差異顯著的腦區數量來確認維度m優化值，如圖3所示。結果發現，m = 1時熵值差異顯著的腦區數量明顯高于m = 2時（P ≤ 0.05），在τ = 2時，差異較明顯，在τ = 3時最為明顯，從τ平均角度上看亦如此。很顯然，m = 1較m = 2更優，在后面的優化過程中將保持m = 1。另外，在τ = 2～3中，r = 0.54～0.6時熵值差異顯著的腦區數量較其他r值要多，即其他2個參數優化將在r = 0.54~0.6、τ = 2～3區間進行。
圖3
以組間熵值差異顯著的腦區數量確定維度m優化值
Figure3.
Optimal value of parameter m estimated by number of brain regions with inter-group significant differences in MSE model
2.1.2 相似度r優化值
保持維度m = 1時，結合ROC曲線和AUC值對相似度r = 0.54～0.6（步長為0.02）進行優化，如圖4所示（τ = 2和τ = 3結果相似，這里只顯示了時間尺度τ = 3的結果）。結果發現，當r = 0.54～0.6時，圖4a左內側額上回和圖4b右頂上回的ROC曲線都處于參考線以上，說明這兩個腦區對致癇側敏感，可作為標志性腦區。然而圖4c左楔葉和圖4d左丘腦的ROC曲線在參考線上下浮動，說明這兩個腦區不具有分類價值，所以不能作為對致癇側敏感的腦區。
圖4
相似系數r對單個腦區分類效果的ROC曲線(m = 1，τ = 3)
a.左內側額上回；b.右頂上回；c.左楔葉；d.左丘腦
Figure4. Effects of similarity coefficient r on classification of a single brain area displayed by ROC curves (m = 1, τ = 3)a. SFGmed.L; b. SPG.R; c. CUN.L; d. THA.L
這個規律也可通過腦區AUC值體現（τ = 2和τ = 3結果相似，這里只給出了時間尺度τ = 3的結果），如表1所示。結果發現，由ROC曲線確認的標志性腦區左內側額上回和右頂上回，其AUC值也較大，且在r = 0.56時值最大，而左楔葉和左丘腦并不存在合適的r值使其AUC值較大。綜上所述，r = 0.56為r的優化值。
表1
不同相似度r下單個腦區側別分類效果的AUC值
Table1.
AUC value of classification of a single brain region by different similarity factor r
2.1.3 時間尺度τ優化值
保持m = 1、r = 0.56下，結合ROC曲線和AUC值對時間尺度 τ = 1～5進行優化，如圖5和表2所示。結果發現，圖5a左腦島和圖5b左中央溝蓋在時間尺度因子 τ = 1～5范圍內，其ROC曲線多數在參考線之上，說明該腦區可作為對致癇側敏感的腦區。然而，圖5c左楔葉和圖5d右背外側額上回的ROC曲線在參考線上下浮動，不具備分類價值，所以不能作為對致癇側敏感的腦區。
圖5
時間尺度因子τ對單個腦區分類效果的ROC曲線（m = 1，r = 0.56）
a.左腦島；b.左中央溝蓋；c.左楔葉；d.右背外側額上回
Figure5. Effects of scale factor τ on classification of a single brain area displayed by ROC curves (m = 1, r = 0.56)a. INS.L; b. ROL.L; c. CUN.L; d. SFGdor.R
表2
不同時間尺度τ下單個腦區側別分類的AUC值
Table2.
AUC value of classification of a single brain region by different scale factor τ
表2給出單個腦區側別分類的AUC值，與ROC曲線規律一致。左腦島和左中央溝蓋在τ = 3時AUC值達到最大，而左楔葉和右背外側額上回則不存在合適的 τ使其AUC值較大。綜上所述，τ = 3為 τ的優化值。
綜上所述，得到TLE患者rfMRI信號多尺度熵模型的優化參數：維度m = 1，相似度r = 0.56，時間尺度 τ = 3。
2.2 特征向量
在優化參數（m = 1，r = 0.56，τ = 3）下，計算TLE患者全腦熵值，通過t檢驗獲得P值（P < 0.01），共選出9個標志性腦區，即左中央溝蓋、右中央溝蓋、左內側額上回、左腦島、右腦島、右頂上回、左楔前葉、右尾狀核、右顳極：顳上回。各腦區ROC曲線及腦區在皮層表面的投影如圖6所示。
圖6
標志性腦區
a.在多尺度熵模型優化參數（即
a. ROC curves of 9 bio-marker brain regions at optimized parameters (
將以上腦區AUC值和P值列入表3，可以看出標志性腦區AUC值較大，并且其P值皆小于0.05。
表3
標志性腦區的AUC值和統計值
Table3.
AUC value and significance (P value) in a bio-marker area
最后，9個標志性腦區在 τ = 1～5時的組間多尺度熵值變化如圖7所示。從圖中看出在 τ = 3時，腦區的組間熵值都表現出最大差異，也驗證了 τ = 3為優化值。
圖7
9個標志性腦區的組間多尺度熵值變化
a.左中央溝蓋；b.右中央溝蓋；c.左內側額上回；d.左腦島；e.右腦島；f.右頂上回；g.左楔前葉；h.右尾狀核；i.右顳極：顳上回
Figure7. Inter-group MSE values change of 9 bio-marker brain regionsa. ROL.L; b. ROL.R; c. SFGmed.L; d. INS.L; e. INS.R; f. SPG.R; g. PCUN.L; h. CAU.R; i. TPOsup.R
2.3 SVM分類效果
將9個標志性腦區的優化熵值作為特征向量，按照3∶2比例將數據隨機分為兩組數據，輸入到SVM進行定側，并用LOOCV進一步驗證，觀察其定側準確度。
經過20次隨機分組計算，定側精度如表4所示，表中精度（
s）為：（92.% ± 7.29）%。留一交叉驗證法的定側結果如表5所示，平均精度為95%。從表4、5中發現，利用優化MSE模型對TLE患者定側精度較為理想。
表4
隨機分組的分類精度（%）
Table4.
Classification rate (CR) of side classification by random validation (%)
表5
LOOCV的分類結果
Table5.
Classification result by LOOCV
3 討論
為進一步體現優化熵模型對區分TLE患者致癇側的有效性，我們在本文中相同TLE患者數據上計算了由Pearson相關系數表達的傳統功能連接（記成FC），由本文采用的敏感性分析（即ROC曲線和AUC值）獲得8個組間差異顯著腦區的FC為特征向量，并結合SVM進行10次隨機交叉驗證，分類結果如表6所示，其分類準確度（
s）為（68.75 ± 12.81）%。與熵模型相比，FC的分類準確度波動范圍較大即穩定性較差，不如本文提出的優化熵模型的分類效果。
表6
基于皮爾森相關性的SVM機器學習分類精度
Table6.
Classification accuracy in machine learning of SVM by FC
表7[9,11-12,37-39]是本文研究方法與其他方法對TLE患者致癇側的分類效果對比。從中可以看出熵模型優化方法有助于分析rfMRI數據，為TLE患者定側研究提供了一種新方法。
表7
本文與其他分析方法定側準確率對比
Table7.
Comparison of classification rate by this work and other studies in localizing epileptogenic hemisphere
目前比較明確的TLE患者病灶包含島葉、雙側扣帶回和楔前葉[12,38-40]。通過分析TLE患者BOLD信號的變化以及對靜息態腦功能網絡的影響，發現楔前葉、左右額葉和頂葉中的BOLD信號顯著減少，聽覺網絡中的腦島及顳上回的功能連接強度降低。此外，還發現TLE患者背側注意網絡中的頂上葉表現出的功能連接強度存在差異[41-42]。利用基于體素的形態測量法分析不同發作側TLE患者的結構像，發現顳葉、額葉和枕葉區域的皮質厚度發生變化。通過容積磁共振成像技術對TLE患者和健康對照組的海馬不對稱性分析得到94%的分類準確度[43]，但TLE患者的結構像變化不局限于海馬體，可能會延伸到顳葉周圍的皮質區域和皮質下結構[44-45]。這些研究結果為本文的標志性腦區提供了支撐。
盡管文中表明優化的多尺度熵模型可較為有效地提取對TLE患者致癇側敏感的功能影像學標記，為優化多尺度熵模型的參數選擇提供了一定的參考性。但無法確定本文所選參數是否適用于不同時間長度的癲癇腦部rfMRI信號，還需要借助其他癲癇腦部rfMRI數據集驗證。文中進行機器學習分類的樣本量較少，雖獲得較為理想的分類準確率，但仍可能存在泛化能力不足的問題。目前無法確定多尺度熵模型結合機器學習在充足樣本中能否取得更好的分類效果，需要在此基礎上進一步搜集數據深入研究。
4 結論
多尺度熵模型應用于TLE患者rfMRI，其首要問題是如何客觀選擇模型的計算參數。為此，以定位致癇側為例，論文提出借助敏感性分析的ROC曲線和AUC值對多尺度熵模型參數進行優化，使其選擇具有一定的客觀性。為檢驗優化效果，用優化參數的熵模型尋求對致癇側敏感的標志性腦區，以其熵值構成特征向量輸入支持向量機進行定側，獲得了較高準確率。分析結果發現9個標志性腦區，即左中央溝蓋、右中央溝蓋、左內側額上回、左腦島、右腦島、右頂上回、左楔前葉、右尾狀核、右顳極：顳上回。將腦區熵值作為特征向量對致癇側不同的TLE患者的分類準確率可達到95%，說明優化多尺度熵模型能夠較為有效地區分TLE患者致癇側。
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
引言
癲癇是由大腦神經元異常放電引起的中樞神經系統功能異常的慢性疾病，嚴重影響著患者生活和生命質量，且病變可在大規模腦網絡中產生分布效應，導致與原發病灶相連的遠端腦區發生代謝功能障礙或者結構改變[1]。手術治療是難治性癲癇的首選方法，70%以上的局灶性癲癇患者手術治療后病情得到緩解，顳葉癲癇（temporal lobe epilepsy，TLE）患者表現出最積極的治療效果[2-3]。因此，癲癇致癇側、致癇灶和癲癇傳播路徑等定位是術前的重要準備，定側是基礎[4]。
臨床上，對于癲癇疾病診斷以及致癇灶定位常采用腦電圖和常規結構磁共振成像（magnetic resonance imaging，MRI）技術。腦電圖是臨床上定位致癇灶和致癇側的金標準，其通過監測始發異常神經放電實現定位。然而，顱內電極或立體定位電極的插入會不可避免地造成腦組織損傷，并且腦電存在空間分辨率低、靈敏度不足、觀察范圍有限和費用昂貴等缺點[5]。此外，常規MRI作為TLE患者術前定位最常用的檢查方法，通過對海馬形態的觀察和海馬體體積測量能夠對部分TLE患者作出定位診斷[6]。海馬陽性是致癇側的典型標識，但只有海馬區神經元丟失達50%機體形態才會表現陽性的異常標記[7]，所以，常規MRI掃描對致癇灶和致癇側定位診斷的檢出率（33%）較低[8]。
近年來，靜息態功能磁共振成像（resting-state functional MRI，rfMRI）以其無損、快捷和廣泛的患者適應性等優勢，成為檢測癲癇引起異常腦功能的先進技術。研究表明，不同致癇側的TLE患者rfMRI存在較大差異[9-10]。Barron等[11]發現不同致癇側的TLE患者丘腦和顳葉之間功能連接強度存在差異；Yang等[12]分析TLE患者的靜息態功能連接，發現該方法選取的特征量能夠較為有效地區分TLE患者的致癇側。rfMRI不僅可以從腦功能水平揭示疾病腦區病理變化，還能增加手術的安全性和準確性，對神經疾病術后及療效監測有較大的臨床應用價值[13-15]。
以往，熵模型廣泛應用于癲癇腦電信號分析。為從神經信號中提取更多有效信息，Costa等[16]提出多尺度熵（multi-scale entropy，MSE）模型將熵和尺度聯合共同表征生理系統在不同狀態下所表現的復雜特性。MSE模型被廣泛用于分析具有多尺度和非線性特點的腦電[17]、fMRI[18-19]、功能性近紅外光譜（functional near-infrared spectroscopy，fNIRS）[20]和心腦信號[21-23]。Yang等[24]通過MSE模型發現年輕人和老年人的血氧水平依賴（blood oxygen level dependent，BOLD）信號復雜度存在差異，Li等[20]通過MSE模型發現輕度認知障礙患者、阿爾茨海默病患者和健康人的fNIRS信號復雜度也存在差異。然而，MSE模型的參數選擇對信號本身存在較大依賴性，以往多依賴經驗選取參數，具有一定的主觀性。對rfMRI信號來說，常采用m = 1、r = 0.5[19]或m = 1、r = 0.35[24]等。
TLE患者rfMRI信號的非線性顯著高于健康人[25-26]，因此，應用MSE模型研究TLE患者的rfMRI復雜度，或可為患者的檢測或治療提供新思路。本文以rfMRI致癇側定位為例，試圖用敏感性分析的接受者操作特征曲線（receiver operating characteristic，ROC）及其曲線下面積AUC值（area under ROC curve，AUC）優化MSE模型的計算參數。為檢驗優化效果，通過優化后的熵模型提取對致癇側敏感的功能像標記，并將其熵值作為特征向量輸入至可適應小樣本的機器學習模型——支持向量機（support vector machine，SVM）定側，以其精度評估優化有效性，試圖解決MSE參數選擇依賴經驗和數據特征、缺少統一客觀依據的計算瓶頸，為熵值應用于TLE患者rfMRI提供客觀依據。最后，由于SVM是有監督學習模型，需要在正確標注的樣本上運行，所以，本文借助結構像海馬陽性提示的致癇側信息進行有監督的機器學習。整體分析框架如圖1所示。
圖1
分析流程圖
Figure1.
Study flowchart
1 數據與方法
1.1 數據采集
本文中20名TLE患者rfMRI數據來自南京大學醫學院附屬金陵醫院醫學影像科，其中男性9例，女性11例，年齡19～33（25.7±3.7）歲，接受常規的術前評估以定位癲癇放電區域。其中結構像海馬陽性左側和右側的患者各10例（男性，左側4例，右側5例；女性，左側6例，右側5例），命名為左側組和右側組。
數據采集使用Siemens 3.0T高場強超導MR掃描儀，經單次激發GE-EPI序列采集BOLD-fMRI數據，TR = 2 000 ms，TE = 30 ms，翻轉角度為90°，矩陣64×64，掃描視野240 mm×240 mm，層厚4.00 mm，層間距0.40 mm[27]。本研究獲得南京軍區總醫院醫學倫理委員會批準，患者知情并且簽署知情同意書。
1.2 數據預處理
使用FMRIB軟件庫中的工具進行rfMRI數據預處理[28-29]，處理流程和質量檢測（頭動控制量 < 0.2 mm）在美國哈佛醫學院Martinos醫學影像中心完成[30-31]。處理平臺CPU為Intel Xeon Sliver 4112 × 16，64核，操作系統Centos7.6，每個被試者的數據預處理時間約為15 h。
預處理步驟：① 將每個被試采集的前4個時間點數據移除，使信號更加穩定；② 使用SPM2校正切片采集的每個腦體積時間偏移（Wellcome Department of Cognitive Neurology，London，英國）；③ 利用FSL軟件（http://fsl.fmrib.ox.ac.uk/fsl）對頭部運動進行剛體平移和旋轉[32]；④ 全腦平均信號的歸一化，將信號配準到蒙特利爾神經研究所（Montreal Neurological Institute，MNI）的標準空間，重采樣到2 mm各向同性體素大小；⑤ 帶通時間濾波（0.01～0.08 Hz）。
預處理后，使用win10操作系統MatlabR2018b版本的DPABI將數據投射到自動解剖標記（the anatomical automatic labeling atlas，AAL1）上，每個被試者的數據處理時間約為20 min，將每個被試全腦分成116個腦區，每個腦區有246個時間點，時間長度為490 s。
1.3 多尺度熵
1.3.1 多尺度熵模型
多尺度熵模型通過評估量化生物醫學信號的“復雜度”，能夠有效分析信號的動態變化，因此在神經信號分析中具有應用價值，特別是癲癇這種動態疾病。
多尺度熵的計算步驟：
1）將采集長度為N的時間序列進行粗粒化，重新構建一個時間序列如式（1）所示。粗粒化（非重疊式）是根據選取的時間尺度大小，對連續的時間點取平均值，重新構建時間序列。
|
式中為時間尺度。
2）定義多尺度熵模型參數m，其中m為整數，代表比較向量的長度。按順序構建一組m維矢量，即， 。
3）定義多尺度熵模型參數公差閾值（即相似系數r），其中r為實數，代表匹配閾值。定義Ui與Uj之間的距離 為兩個矢量中元素之間距離的最大值，即 ，其中i ≠ j，統計當 時 小于r時的數目b，并將b與距離總數 的比值，記做 ，如式（2）所示。
|
4）將維數m加1，在重新設定維度m的情況下重復步驟2、3，將會得到。
5）理論上，多尺度熵計算如式（3）所示：
|
當N有限時，式（3）的計算模型如式（4）所示：
|
綜上所述，多尺度熵中的3個參數分別為維度m、相似度r和時間尺度τ。若計算所得熵值越大，表明復雜度越大，反之亦然。
1.3.2 熵模型參數優化準側
多尺度熵是在不同時間尺度值下計算時間序列的樣本熵，所以，多尺度熵模型中維度m和相似度r的選取與樣本熵參數選取相似。因為多尺度熵的3個參數之間相互影響，則在單個參數優化過程中，其余兩個參數也須參與優化過程，即3個參數同時參與一個參數的優化過程。本文中，首先利用組間熵值差異較大（P ≤ 0.05）的腦區數量確定維度m和另外兩個參數即將優化的區間，然后，結合ROC曲線和AUC值依次對相似度r和時間尺度 τ 進行優化。多尺度熵模型的參數優化流程見圖2。
圖2
模型參數優化流程圖
Figure2.
Flowchart of parameters optimization in a multi-scale entropy model
原則上，熵評估的準確度會隨著維度m和m + 1的矢量（時間序列）匹配數量的增加而提高，選擇較小的維度m和較大的相似度r也能夠增加矢量匹配數量提高準確度。基于時間長度10m～20m（m參數維度）的BOLD信號可較為準確地計算熵值[33]。根據本文中數據的時間點長度可以選取m = 1或2。
根據經驗，相似度r范圍可選在0.05～0.6。但對文中數據分析發現：在r = 0.05～0.3時，腦區熵值會出現無效值，這是由于選取相似度r較小，計算過程中不存在與之相匹配的矢量。為保證分析精確性，相似度r選在0.3～0.6（步長0.02）范圍內，時間尺度也不宜選擇過大，本文選取范圍是 = 1～5。
本文用敏感性分析的ROC曲線（反映分類效果敏感性和特異性連續變量的綜合指標）和AUC值（ROC曲線下方的面積大小）來評估MSE模型參數的優化，使參數有利于分析癲癇rfMRI信號特征，使得參數選擇具有一定的統計特性和客觀性。
1.3.3 對發作側敏感的腦區
利用MSE模型計算全腦（AAL1）中90個大腦區域的rfMRI信號得到腦區熵值，并通過t檢驗對比兩組患者腦區熵值的差異性，獲得腦區P值。在P ≤ 0.05范圍內，選取腦區將其作為對致癇側敏感的功能影像標記，并利用ROC曲線和AUC值（IBM SPSS Statistics21，美國）對標記進行驗證。最后用BrainNet Viewer[34]（http://www.nitrc.org/projects/bnv/）將標記進行可視化。
1.4 機器學習
1.4.1 特征向量
在優化參數下，將對致癇側敏感的腦區熵值作為特征向量，輸入到支持向量機中。
1.4.2 支持向量機
由Vapnik提出的SVM是一種有監督的機器學習方法，以訓練誤差作為優化問題的約束條件，并以置信范圍最小化為優化目標，是一種基于結構風險最小化準則的學習方法，不易受數據維數和樣本數量的限制，這為其在小樣本中的應用提供了依據[35-36]。
文中將20名患者的腦區多尺度熵值作為特征向量，按照3∶2比例分成12名訓練集和8名測試集，作為機器學習的基本輸入信息進行數據歸一化。選擇徑向基函數（radial basis function，RBF）做為核函數，選取比例參數g = 代替核參數，形成一組參數對（C，g）。其中，參數C和g的取值范圍設定為［? 10，10］（步長0.2），并用網格搜索方法（交叉驗證法）尋找誤差最小的懲罰系數C和核參數為最優參數，創建訓練模型，對測試集數據進行計算。
1.4.3 留一交叉驗證
留一交叉驗證（leave one out cross-validation，LOOCV）常被用來檢測機器學習模型的準確性，其樣本利用率較高，能夠從有限的樣本中盡可能地獲取更多有效信息，所以，適用于小樣本分析。如果有N個樣本，其中N ? 1個為訓練樣本，剩余1個為測試樣本，用N次結果的平均分類精度來評估模型的準確性（文中，N = 20）。
2 結果
2.1 MSE模型的優化參數
調節3個參數（維度m、相似度r和時間尺度τ）計算腦區熵值，通過t檢驗分析左、右致癇側TLE患者間的腦區熵值，從而獲得具有統計學意義（P ≤ 0.05）的腦區和優化參數范圍，并利用ROC曲線和AUC值對其驗證。
2.1.1 維度m優化值
以時間尺度τ = 1～5（步長1）生成圖，每張圖以相似度r = 0.3～0.6（步長為0.02）為變量，對比維度m = 1和m = 2時組間熵值差異顯著的腦區數量來確認維度m優化值，如圖3所示。結果發現，m = 1時熵值差異顯著的腦區數量明顯高于m = 2時（P ≤ 0.05），在τ = 2時，差異較明顯，在τ = 3時最為明顯，從τ平均角度上看亦如此。很顯然，m = 1較m = 2更優，在后面的優化過程中將保持m = 1。另外，在τ = 2～3中，r = 0.54～0.6時熵值差異顯著的腦區數量較其他r值要多，即其他2個參數優化將在r = 0.54~0.6、τ = 2～3區間進行。
圖3
以組間熵值差異顯著的腦區數量確定維度m優化值
Figure3.
Optimal value of parameter m estimated by number of brain regions with inter-group significant differences in MSE model
2.1.2 相似度r優化值
保持維度m = 1時，結合ROC曲線和AUC值對相似度r = 0.54～0.6（步長為0.02）進行優化，如圖4所示（τ = 2和τ = 3結果相似，這里只顯示了時間尺度τ = 3的結果）。結果發現，當r = 0.54～0.6時，圖4a左內側額上回和圖4b右頂上回的ROC曲線都處于參考線以上，說明這兩個腦區對致癇側敏感，可作為標志性腦區。然而圖4c左楔葉和圖4d左丘腦的ROC曲線在參考線上下浮動，說明這兩個腦區不具有分類價值，所以不能作為對致癇側敏感的腦區。
圖4
相似系數r對單個腦區分類效果的ROC曲線(m = 1，τ = 3)
a.左內側額上回；b.右頂上回；c.左楔葉；d.左丘腦
Figure4. Effects of similarity coefficient r on classification of a single brain area displayed by ROC curves (m = 1, τ = 3)a. SFGmed.L; b. SPG.R; c. CUN.L; d. THA.L
這個規律也可通過腦區AUC值體現（τ = 2和τ = 3結果相似，這里只給出了時間尺度τ = 3的結果），如表1所示。結果發現，由ROC曲線確認的標志性腦區左內側額上回和右頂上回，其AUC值也較大，且在r = 0.56時值最大，而左楔葉和左丘腦并不存在合適的r值使其AUC值較大。綜上所述，r = 0.56為r的優化值。
表1
不同相似度r下單個腦區側別分類效果的AUC值
Table1.
AUC value of classification of a single brain region by different similarity factor r
2.1.3 時間尺度τ優化值
保持m = 1、r = 0.56下，結合ROC曲線和AUC值對時間尺度 τ = 1～5進行優化，如圖5和表2所示。結果發現，圖5a左腦島和圖5b左中央溝蓋在時間尺度因子 τ = 1～5范圍內，其ROC曲線多數在參考線之上，說明該腦區可作為對致癇側敏感的腦區。然而，圖5c左楔葉和圖5d右背外側額上回的ROC曲線在參考線上下浮動，不具備分類價值，所以不能作為對致癇側敏感的腦區。
圖5
時間尺度因子τ對單個腦區分類效果的ROC曲線（m = 1，r = 0.56）
a.左腦島；b.左中央溝蓋；c.左楔葉；d.右背外側額上回
Figure5. Effects of scale factor τ on classification of a single brain area displayed by ROC curves (m = 1, r = 0.56)a. INS.L; b. ROL.L; c. CUN.L; d. SFGdor.R
表2
不同時間尺度τ下單個腦區側別分類的AUC值
Table2.
AUC value of classification of a single brain region by different scale factor τ
表2給出單個腦區側別分類的AUC值，與ROC曲線規律一致。左腦島和左中央溝蓋在τ = 3時AUC值達到最大，而左楔葉和右背外側額上回則不存在合適的 τ使其AUC值較大。綜上所述，τ = 3為 τ的優化值。
綜上所述，得到TLE患者rfMRI信號多尺度熵模型的優化參數：維度m = 1，相似度r = 0.56，時間尺度 τ = 3。
2.2 特征向量
在優化參數（m = 1，r = 0.56，τ = 3）下，計算TLE患者全腦熵值，通過t檢驗獲得P值（P < 0.01），共選出9個標志性腦區，即左中央溝蓋、右中央溝蓋、左內側額上回、左腦島、右腦島、右頂上回、左楔前葉、右尾狀核、右顳極：顳上回。各腦區ROC曲線及腦區在皮層表面的投影如圖6所示。
圖6
標志性腦區
a.在多尺度熵模型優化參數（即
a. ROC curves of 9 bio-marker brain regions at optimized parameters (
將以上腦區AUC值和P值列入表3，可以看出標志性腦區AUC值較大，并且其P值皆小于0.05。
表3
標志性腦區的AUC值和統計值
Table3.
AUC value and significance (P value) in a bio-marker area
最后，9個標志性腦區在 τ = 1～5時的組間多尺度熵值變化如圖7所示。從圖中看出在 τ = 3時，腦區的組間熵值都表現出最大差異，也驗證了 τ = 3為優化值。
圖7
9個標志性腦區的組間多尺度熵值變化
a.左中央溝蓋；b.右中央溝蓋；c.左內側額上回；d.左腦島；e.右腦島；f.右頂上回；g.左楔前葉；h.右尾狀核；i.右顳極：顳上回
Figure7. Inter-group MSE values change of 9 bio-marker brain regionsa. ROL.L; b. ROL.R; c. SFGmed.L; d. INS.L; e. INS.R; f. SPG.R; g. PCUN.L; h. CAU.R; i. TPOsup.R
2.3 SVM分類效果
將9個標志性腦區的優化熵值作為特征向量，按照3∶2比例將數據隨機分為兩組數據，輸入到SVM進行定側，并用LOOCV進一步驗證，觀察其定側準確度。
經過20次隨機分組計算，定側精度如表4所示，表中精度（s）為：（92.% ± 7.29）%。留一交叉驗證法的定側結果如表5所示，平均精度為95%。從表4、5中發現，利用優化MSE模型對TLE患者定側精度較為理想。
表4
隨機分組的分類精度（%）
Table4.
Classification rate (CR) of side classification by random validation (%)
表5
LOOCV的分類結果
Table5.
Classification result by LOOCV
3 討論
為進一步體現優化熵模型對區分TLE患者致癇側的有效性，我們在本文中相同TLE患者數據上計算了由Pearson相關系數表達的傳統功能連接（記成FC），由本文采用的敏感性分析（即ROC曲線和AUC值）獲得8個組間差異顯著腦區的FC為特征向量，并結合SVM進行10次隨機交叉驗證，分類結果如表6所示，其分類準確度（s）為（68.75 ± 12.81）%。與熵模型相比，FC的分類準確度波動范圍較大即穩定性較差，不如本文提出的優化熵模型的分類效果。
表6
基于皮爾森相關性的SVM機器學習分類精度
Table6.
Classification accuracy in machine learning of SVM by FC
表7[9,11-12,37-39]是本文研究方法與其他方法對TLE患者致癇側的分類效果對比。從中可以看出熵模型優化方法有助于分析rfMRI數據，為TLE患者定側研究提供了一種新方法。
表7
本文與其他分析方法定側準確率對比
Table7.
Comparison of classification rate by this work and other studies in localizing epileptogenic hemisphere
目前比較明確的TLE患者病灶包含島葉、雙側扣帶回和楔前葉[12,38-40]。通過分析TLE患者BOLD信號的變化以及對靜息態腦功能網絡的影響，發現楔前葉、左右額葉和頂葉中的BOLD信號顯著減少，聽覺網絡中的腦島及顳上回的功能連接強度降低。此外，還發現TLE患者背側注意網絡中的頂上葉表現出的功能連接強度存在差異[41-42]。利用基于體素的形態測量法分析不同發作側TLE患者的結構像，發現顳葉、額葉和枕葉區域的皮質厚度發生變化。通過容積磁共振成像技術對TLE患者和健康對照組的海馬不對稱性分析得到94%的分類準確度[43]，但TLE患者的結構像變化不局限于海馬體，可能會延伸到顳葉周圍的皮質區域和皮質下結構[44-45]。這些研究結果為本文的標志性腦區提供了支撐。
盡管文中表明優化的多尺度熵模型可較為有效地提取對TLE患者致癇側敏感的功能影像學標記，為優化多尺度熵模型的參數選擇提供了一定的參考性。但無法確定本文所選參數是否適用于不同時間長度的癲癇腦部rfMRI信號，還需要借助其他癲癇腦部rfMRI數據集驗證。文中進行機器學習分類的樣本量較少，雖獲得較為理想的分類準確率，但仍可能存在泛化能力不足的問題。目前無法確定多尺度熵模型結合機器學習在充足樣本中能否取得更好的分類效果，需要在此基礎上進一步搜集數據深入研究。
4 結論
多尺度熵模型應用于TLE患者rfMRI，其首要問題是如何客觀選擇模型的計算參數。為此，以定位致癇側為例，論文提出借助敏感性分析的ROC曲線和AUC值對多尺度熵模型參數進行優化，使其選擇具有一定的客觀性。為檢驗優化效果，用優化參數的熵模型尋求對致癇側敏感的標志性腦區，以其熵值構成特征向量輸入支持向量機進行定側，獲得了較高準確率。分析結果發現9個標志性腦區，即左中央溝蓋、右中央溝蓋、左內側額上回、左腦島、右腦島、右頂上回、左楔前葉、右尾狀核、右顳極：顳上回。將腦區熵值作為特征向量對致癇側不同的TLE患者的分類準確率可達到95%，說明優化多尺度熵模型能夠較為有效地區分TLE患者致癇側。
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
