基于功率譜及有限穿越可視圖的癲癇腦電信號分析算法

王若凡，劉 靜，王 江，于海濤，曹亦賓

(1.天津職業(yè)技術(shù)師范大學(xué) 信息技術(shù)工程學(xué)院，天津 300222；2.唐山市工人醫(yī)院 神經(jīng)內(nèi)科，河北 唐山 064300； 3.天津大學(xué) 電氣與自動化工程學(xué)院，天津 300072)

(*通信作者電子郵箱m18522619839@163.com)

基于功率譜及有限穿越可視圖的癲癇腦電信號分析算法

王若凡1，劉 靜2，王 江3，于海濤3，曹亦賓2*

(1.天津職業(yè)技術(shù)師范大學(xué) 信息技術(shù)工程學(xué)院，天津 300222；2.唐山市工人醫(yī)院 神經(jīng)內(nèi)科，河北 唐山 064300； 3.天津大學(xué) 電氣與自動化工程學(xué)院，天津 300072)

(*通信作者電子郵箱m18522619839@163.com)

針對可視圖(VG)算法存在噪聲魯棒性差的問題，提出一種改進(jìn)的有限穿越可視圖(LPVG)建網(wǎng)方法。該算法基于可視圖(VG)算法的可視性準(zhǔn)則，并設(shè)定有限穿越視距，將時間序列中滿足條件的點(diǎn)連接起來，從而將時間序列映射為網(wǎng)絡(luò)。首先，對LPVG算法進(jìn)行性能分析；然后，將LPVG算法結(jié)合功率譜密度(PSD)算法應(yīng)用到癲癇發(fā)作前、中、后腦電信號的識別上；最后，提取三種狀態(tài)下癲癇腦電信號的LPVG網(wǎng)絡(luò)特征參數(shù)，研究癲癇對網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)的影響。仿真結(jié)果表明，與VG和水平穿越可視圖(HVG)相比，雖然LPVG算法的時間復(fù)雜度較高，但是LPVG對信號中的噪聲具有較強(qiáng)的魯棒性：分別對周期、隨機(jī)、分形和混沌四種時間序列進(jìn)行LPVG建網(wǎng)，發(fā)現(xiàn)隨著噪聲強(qiáng)度增大，LPVG網(wǎng)絡(luò)聚類系數(shù)的波動率均為最低，分別為6.73%、0.05%、0.99%和3.20%。接下來對腦電信號的PSD和LPVG建網(wǎng)分析結(jié)果表明，癲癇發(fā)作中，PSD值在delta頻帶下顯著增強(qiáng)，而在theta頻帶下顯著降低；LPVG網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)有所改變，網(wǎng)絡(luò)中各模塊的獨(dú)立性有所提高，網(wǎng)絡(luò)的平均路徑長度增大，復(fù)雜度降低。所提的功率譜密度和有限穿越可視圖算法能夠有效表征癲癇前、中、后三種狀態(tài)下的腦電信號能量分布和單通道信號可視化后的網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)的異常，為癲癇的病理研究和臨床診斷提供幫助。

腦電信號；癲癇；功率譜密度；有限穿越可視圖；復(fù)雜網(wǎng)絡(luò)

0 引言

癲癇是一種嚴(yán)重危害人類安全健康的常見慢性腦部疾病，全世界大約0.5%的人患有癲癇，使運(yùn)動感覺、自主神經(jīng)、意識和精神狀態(tài)等反復(fù)出現(xiàn)不同程度的障礙。癲癇發(fā)作時表現(xiàn)為大腦神經(jīng)元興性增高以及過度同步化放電，并導(dǎo)致短暫性中樞神經(jīng)系統(tǒng)功能失常[1-2]，患者承受了巨大的痛苦和危險，極大影響了患者及其家庭的生活質(zhì)量。

腦電圖(ElectroEncephaloGram, EEG)是在頭皮上通過電極記錄的腦細(xì)胞群自發(fā)性、節(jié)律性的電活動軌跡[3]。當(dāng)大腦出現(xiàn)病理性或功能性改變時，EEG就會發(fā)生相應(yīng)變化[4-5]。此外，腦電圖具有時間分辨率高、價格低廉、獲取容易且可進(jìn)行實(shí)時監(jiān)測等優(yōu)點(diǎn)，日益成為癲癇臨床診斷的重要輔助工具，在臨床醫(yī)學(xué)和腦科學(xué)的研究中都起著不可替代的重要作用[5-7]。癲癇診斷目前大都由醫(yī)生通過視覺檢測患者的腦電圖，根據(jù)經(jīng)驗(yàn)診斷完成的。人工視覺檢測存在費(fèi)時、效率低、且缺乏標(biāo)準(zhǔn)的制約等不利因素，因此，癲癇腦電的自動識別與診斷研究具有重大的臨床應(yīng)用價值。

針對這一問題，眾多學(xué)者致力于癲癇EEG的自動檢測技術(shù)領(lǐng)域的研究，并提出很多自動檢測方法。其中最常見的分析方法主要有時域分析、頻域分析、復(fù)雜度分析等[8-16]。這些自動方法診斷主要是對腦電信號進(jìn)行分析，提取癲癇疾病的特征向量，以實(shí)現(xiàn)對健康人群與癲癇患者的區(qū)分以及癲癇不同時期的區(qū)分。大部分的參考文獻(xiàn)基于對不同的EEG通道進(jìn)行譜分析[8-10]。Mporas等[11]通過使用6階自回歸(AutoRegression, AR)模型估計癲癇EEG的功率譜密度并以此為特征進(jìn)行分類分析，準(zhǔn)確率高達(dá)近90%。此外，樣本熵、多尺度熵、模糊熵、相空間分析等非線性混沌算法也被用來提取腦電特征，用于癲癇的自動診斷研究[12-16]。

近年來，基于復(fù)雜網(wǎng)絡(luò)的非線性時間序列的動力學(xué)分析受到了廣泛關(guān)注[17-20]。如Zhang等[20]在相空間內(nèi)將準(zhǔn)周期時間序列的每一周期作為網(wǎng)絡(luò)的一個節(jié)點(diǎn)，每兩個準(zhǔn)周期間的距離與所選取的閾值共同決定節(jié)點(diǎn)之間是否存在連邊，以此將心電信號轉(zhuǎn)化為復(fù)雜網(wǎng)絡(luò)用于心電信號的分析。Gao等[21]將一維時間序列重構(gòu)到高維的相空間中，由此構(gòu)成對應(yīng)原時間序列的復(fù)雜網(wǎng)絡(luò)，并將此相空間重構(gòu)方法應(yīng)用于識別與分析多相流的流形。Wang等[22]通過相干性分析將多通道Alzheimer癥腦電映射為復(fù)雜腦網(wǎng)絡(luò)，通過網(wǎng)絡(luò)特征探究Alzheimer腦功能的異常。

2008年Lacasa等[23]提出了一種全新的一維時間序列映射為復(fù)雜網(wǎng)絡(luò)的方法：可視圖(Visibility Graph, VG)算法。其思想是將離散時間序列中的數(shù)據(jù)點(diǎn)作為網(wǎng)絡(luò)節(jié)點(diǎn)，數(shù)據(jù)點(diǎn)之間滿足可視化準(zhǔn)則的連線定義為網(wǎng)絡(luò)連接?？梢晥D算法已成功應(yīng)用于大腦神經(jīng)障礙分析、心率分析、飛機(jī)湍射流分析、海洋潮汐記錄分析、金融經(jīng)濟(jì)分析[24-29]。隨后，Luque等[30]提出了水平可視圖(Horizontal Visibility Graph, HVG)算法?？梢晥D與水平可視圖建網(wǎng)方法具有算法簡單易實(shí)現(xiàn)、計算效率高等優(yōu)點(diǎn)。但是這兩種算法對不同類別的時間序列的識別能力尚存在問題，且對噪聲的魯棒性差?；诖?，周婷婷等[31]在此基礎(chǔ)上提出了有限穿越可視圖(Limited Penetrable Visibility Graph， LPVG)算法，并將其應(yīng)用于油氣水三相流電導(dǎo)波動序列的識別中。Pei等[32]采用LPVG算法成功地區(qū)分開了針刺前、中、后的腦電信號。本文采用功率譜密度分析和有限穿越可視圖(LPVG)算法對癲癇發(fā)作前、中、后階段的腦電信號進(jìn)行分析，旨在獲取單通道的癲癇腦電序列的譜特征以及EEG序列映射成復(fù)雜網(wǎng)絡(luò)的網(wǎng)絡(luò)特征，用于癲癇腦電信號的自動識別與診斷，為癲癇的病理研究和臨床診斷提供參考。

1 癲癇腦電的采集與預(yù)處理

研究所需癲癇病例由河北省唐山市工人醫(yī)院神經(jīng)內(nèi)科提供。目前，臨床實(shí)驗(yàn)獲取癲癇病例已達(dá)到70余人。本研究在已獲取的病例群中篩選特性相同或相近的5例中度異常病例進(jìn)行算法研究。5例病例均為EP大發(fā)作，具體癥狀表現(xiàn)為右側(cè)頂中后顳持續(xù)大量中高幅尖波尖慢。所有研究病例均為自愿參加，且已簽署知情協(xié)議。

臨床采集腦電實(shí)驗(yàn)中，電極帽電極位置按照國際導(dǎo)聯(lián)10-20系統(tǒng)安放，其中，耳垂作為電極的參考地。實(shí)驗(yàn)設(shè)備的采樣頻率為256 Hz，硬件濾波器的為0.5～100 Hz。

實(shí)驗(yàn)過程中，對實(shí)驗(yàn)病例連續(xù)采集記錄腦電信號10 min。從每位實(shí)驗(yàn)病例的每通道EEG信號中截取發(fā)作前(R1)、中(R2)、后(R3)這三段長度為8 s的數(shù)據(jù)(數(shù)據(jù)長度為256×8=2 048個采樣點(diǎn))進(jìn)行分析。采用0～30 Hz的帶通濾波器分別對每段數(shù)據(jù)進(jìn)行預(yù)處理，并采用小波變換提取EEG信號4個子頻帶：delta(0～4 Hz)、theta(4～8 Hz)、alpha(8～15 Hz)以及beta(15～30 Hz)。圖1為theta頻帶下發(fā)作前、中、后三種狀態(tài)下的腦電信號示意圖，可以發(fā)現(xiàn)相比發(fā)作前、后，癲癇發(fā)作中的腦電信號出現(xiàn)異常。

圖1 theta頻帶下三種狀態(tài)的腦電信號

2 分析方法

本文將采用功率譜密度分析和有限穿越可視圖算法分別從能量和網(wǎng)絡(luò)特征角度對癲癇發(fā)作前、中、后三個狀態(tài)下的單通道腦電信號進(jìn)行分析，具體過程如圖2所示。

2.1 功率譜密度算法

本文采用線性自回歸(AR)模型估計法來估計EEG信號的功率譜密度。

AR模型如式(1)所示：

(1)

其中：x(n)為待分析時間序列，u(n)為外加白噪聲序列，p是AR模型的階數(shù)，本文中取p=10，ak為模型參數(shù)，由Burg算法估計得到，其基本原理是使前向預(yù)測與后向預(yù)測的均方誤差之和最小。

圖2 癲癇腦電分析流程

因此，估計出AR模型參數(shù)后，可以得到序列x(n)的功率譜密度，如式(2)所示：

(2)

2.2 有限穿越可視圖算法

可視圖(VG)算法就是通過可視圖準(zhǔn)則將一維時間序列映射為復(fù)雜網(wǎng)絡(luò)。其思想是將離散時間序列中的數(shù)據(jù)點(diǎn)作為網(wǎng)絡(luò)節(jié)點(diǎn)，數(shù)據(jù)點(diǎn)之間滿足可視化準(zhǔn)則的連線定義為網(wǎng)絡(luò)連接[23]。

可視化準(zhǔn)則如下：長度為N的時間序列{Xi|1,2,…,N}，其中任意兩點(diǎn)Xa與Xc之間存在連接等價于這兩點(diǎn)間任意一點(diǎn)Xb，滿足式(3)：

Xb≤Xc+(Xa-Xc)·(c-b)/(c-a)

(3)

圖3(a)中用黑色垂直的直方條表示時間序列的數(shù)據(jù)點(diǎn)，同時作為網(wǎng)絡(luò)的節(jié)點(diǎn)，直方條的高度表示對應(yīng)的幅值，根據(jù)上述可視化準(zhǔn)則將滿足條件的點(diǎn)連接起來，用實(shí)線表示，構(gòu)成可視圖網(wǎng)絡(luò)。VG網(wǎng)絡(luò)具有如下幾點(diǎn)性質(zhì)：每一個點(diǎn)至少和其前后兩個鄰點(diǎn)相連接；網(wǎng)絡(luò)連接具有無向性；橫軸和縱軸坐標(biāo)尺度變化或者經(jīng)過仿射變換后，可視性保持不變。此外，VG算法能夠繼承原始時間序列的非線性動力學(xué)特征，即應(yīng)用VG算法可將周期時間序列轉(zhuǎn)化成規(guī)則網(wǎng)絡(luò)，隨機(jī)時間序列轉(zhuǎn)化成隨機(jī)網(wǎng)絡(luò)，分形時間序列轉(zhuǎn)化成標(biāo)度網(wǎng)絡(luò)。

水平可視圖(HVG)算法如圖3(b)所示。網(wǎng)絡(luò)的節(jié)點(diǎn)定義與網(wǎng)絡(luò)的連邊方式與VG相同，不同的是，網(wǎng)絡(luò)節(jié)點(diǎn)間的連邊必須是水平直線，且不能穿越其他的直方條。因此，對于同一列時間序列，水平可視圖是可視圖的一個子圖，且節(jié)點(diǎn)的平均度值較小[30]。

有限穿越可視圖(LPVG)算法如圖3(c)所示。定義有限穿越視距為N，當(dāng)且僅當(dāng)兩個節(jié)點(diǎn)所對應(yīng)的直方條之間的連線穿越其他直方條次數(shù)滿足n≤N時，兩節(jié)點(diǎn)之間存在連邊[31]。圖3(c)中實(shí)線表示基于VG算法的連邊，灰色虛線表示新增有限穿越連邊。本文仿真中，有限穿越視距設(shè)為N=5。LPVG算法同樣具有連接無向性與仿射變換或者坐標(biāo)軸尺度變換后可視性不變的特點(diǎn)，且LPVG網(wǎng)絡(luò)具有更強(qiáng)的連接性，長程連接能夠反映序列波動趨勢。因此，LPVG算法對噪聲和干擾具有較強(qiáng)魯棒性，在含噪情況下仍能較好地繼承信號的動力學(xué)特性。

圖3 三種可視圖算法原理示意圖

2.3 網(wǎng)絡(luò)參數(shù)

復(fù)雜網(wǎng)絡(luò)的拓?fù)浣Y(jié)構(gòu)可由以下網(wǎng)絡(luò)參數(shù)[35-37]來表征。

1)節(jié)點(diǎn)度和度分布。

節(jié)點(diǎn)度衡量單個節(jié)點(diǎn)作為功能網(wǎng)絡(luò)中多個連接交叉點(diǎn)的重要程度。網(wǎng)絡(luò)中單個節(jié)點(diǎn)的度值等于該節(jié)點(diǎn)的連接總數(shù)。網(wǎng)絡(luò)的平均度定義為所有節(jié)點(diǎn)度的平均值：

(4)

其中：N是網(wǎng)絡(luò)的節(jié)點(diǎn)個數(shù)，Ki是節(jié)點(diǎn)i的度值，節(jié)點(diǎn)度K的分布情況可以用分布函數(shù)P(K)表示。

2)平均路徑長度。

平均路徑長度是網(wǎng)絡(luò)中所有點(diǎn)對的距離的平均值，度量了網(wǎng)絡(luò)的整體連通性與不同節(jié)點(diǎn)間的集合程度。平均路徑長度定義為某一節(jié)點(diǎn)到任一其他節(jié)點(diǎn)的平均距離：

(5)

其中：dij為節(jié)點(diǎn)i和j間的最短距離，N為網(wǎng)絡(luò)中節(jié)點(diǎn)總個數(shù)。另一個重要的網(wǎng)絡(luò)特征參數(shù)：聚類系數(shù)衡量一個節(jié)點(diǎn)最近鄰節(jié)點(diǎn)之間的連接情況，與平均路徑長度呈現(xiàn)相反趨勢。

3)圖指數(shù)復(fù)雜度。

圖指數(shù)復(fù)雜度(Graph Index Complexity, GIC)用于衡量功能網(wǎng)絡(luò)的復(fù)雜度。定義式如下：

Cλmax=4m(1-m)

(6)

m=[λmax-2 cos(π/(N+1))]/[n-1-2 cos(π/(N+1))]

(7)

其中λmax為功能網(wǎng)絡(luò)鄰接矩陣的最大特征值。在無向二值網(wǎng)絡(luò)中，總有2 cos(π/(n+1))≤λmax≤n-1成立，因此，圖指數(shù)復(fù)雜度在[0,1]區(qū)間變化。Cλmax越大表明網(wǎng)絡(luò)結(jié)構(gòu)越復(fù)雜。

3 模型驗(yàn)證

在本章中，用正弦周期序列、隨機(jī)序列、H=0.5(H為Hurst指數(shù))的布朗分形序列以及Lorenz混沌序列來研究有限穿越可視圖(LPVG)、可視圖(VG)以及水平可視圖(HVG)這三種時間序列建網(wǎng)的可視化算法的性能。此外，對四種序列分別添加10 dB、20 dB、30 dB的高斯白噪聲信號。簡明起見，四種典型時間序列以及添加10 dB高斯白噪聲的時間序列如圖4所示，而其他噪聲強(qiáng)度的序列在圖4中不予表示。

分別采用VG、HVG以及LPVG 3種可視圖算法對四種典型序列進(jìn)行建網(wǎng)，可以發(fā)現(xiàn)3種可視圖算法都具有較好的還原序列的動力學(xué)特性。提取3種可視圖網(wǎng)絡(luò)的聚類系數(shù)以及平均路徑長度，如表1～2所示。可以發(fā)現(xiàn)，隨著噪聲強(qiáng)度的增大，與VG和HVG相比，LPVG算法提取的聚類系數(shù)和平均路徑長度，其波動率均為最小，可見LPVG對含噪信號具有較好的適應(yīng)性，具備更佳的抗干擾能力。

圖4 4種典型仿真時間序列示意圖

表1 3種可視圖算法建網(wǎng)的網(wǎng)絡(luò)聚類系數(shù)對比

Tab.1 Comparison of clustering coefficients of network constructed by three VG algorithms

序列算法網(wǎng)絡(luò)聚類系數(shù)無噪聲10dB20dB30dB波動率/%周期隨機(jī)分形混沌VG0．61740．75420．70070．526316．80HVG0．48870．64570．63540．588327．51LPVG0．60610．64750．62730．54636．73VG0．75210．75580．75130．75160．22HVG0．64530．64670．64520．64520．08LPVG0．65760．65790．65790．65800．05VG0．69380．72810．69890．69622．01HVG0．60390．62580．60580．60321．35LPVG0．62870．64230．63130．62620．99VG0．65980．55900．59260．64019．48HVG0．52240．56440．53600．52733．86LPVG0．69650．64860．68000．69413．20

本文對3種可視圖算法進(jìn)行了時間復(fù)雜度和空間復(fù)雜度分析，計算結(jié)果如表3所示。時間復(fù)雜度是度量算法執(zhí)行的時間長短；而空間復(fù)雜度是度量算法所需存儲空間的大小。

由表3可知，HVG、VG以及LPVG這3種可視圖算法的時間復(fù)雜度依次增大，表明算法的復(fù)雜程度逐漸增大，但計算所需的存儲空間均為O(1)量級。此外，3種可視圖算法的計算時間如圖5(a)所示。仿真時間序列長度設(shè)置為2 500個點(diǎn)。圖5(b)為采用VG與HVG算法對仿真序列的第2 000～2 500點(diǎn)建網(wǎng)時所耗時間的放大。

表2 3種可視圖算法建網(wǎng)的網(wǎng)絡(luò)平均路徑長度對比

表3 3種可視圖算法建網(wǎng)的時間復(fù)雜度

圖5 3種可視圖算法的計算時間比較

由圖5可知，與VG和HVG相比，LPVG算法時間復(fù)雜度遠(yuǎn)大于其他兩種可視圖算法，且仿真消耗時間也遠(yuǎn)遠(yuǎn)大于兩種可視圖算法。

雖然有限穿越可視圖(LPVG)算法的時間復(fù)雜度較高，計算時間較長，且在信號識別能力上存在一定的局限性，但是抗噪能力強(qiáng)。在一定程度上，有限穿越可視圖以犧牲相鄰節(jié)點(diǎn)的相關(guān)性，來抵抗噪聲點(diǎn)的干擾，達(dá)到防止重要節(jié)點(diǎn)被噪聲隔斷的目的，能夠更好地反映時間序列的動力學(xué)特性和波動趨勢。如果有限穿越視距N取值小，局部相關(guān)性并未全部抹殺。有限穿越可視圖在犧牲微觀相關(guān)性的情況下更能保留信號的宏觀信息。

4 實(shí)驗(yàn)結(jié)果分析

本章將第2章提出的算法應(yīng)用到癲癇腦電信號發(fā)作前、中、后狀態(tài)的識別上。由于癲癇腦電信號在臨床采集過程中含有大量噪聲，因此，本章采用功率譜密度(Power Spectrum Density, PSD)以及抗噪能力強(qiáng)的有限穿越可視圖(LPVG)算法對腦電信號(EEG)進(jìn)行分析。

首先將每通道EEG提取為發(fā)作前(R1)、發(fā)作中(R2)、發(fā)作后(R3)三段數(shù)據(jù)。隨后，采用小波分解將三段EEG分解為delta、alpha、theta、beta四個頻帶，并在每個頻帶下采用AR Burg參數(shù)模型法估計分別19通道腦電信號的功率譜密度。最后，對5名患者的PSD計算結(jié)果求平均值，計算結(jié)果如圖6所示。

圖6 癲癇發(fā)作前、中、后的平均功率譜密度地形圖

從圖6中可知，在delta頻帶下(圖6(a))，大腦普遍EEG能量較高，并且在癲癇發(fā)作前與發(fā)作后，腦區(qū)的能量分布有一定的層次性，而在發(fā)作期間，整個腦區(qū)能量異常增強(qiáng)，各腦區(qū)的能量值顯著升高；而在theta頻帶下(圖6(b))，發(fā)作期間整個腦區(qū)的能量異常減少，腦區(qū)之間能量分布的區(qū)域性差異顯著降低，反映出大腦腦區(qū)活動趨于一致性；在alpha頻帶下(圖6(c))，發(fā)作前、中、后三個階段，腦區(qū)能量的變化不明顯；在beta頻帶下(圖6(d))，腦區(qū)能量的分布值很低，且不存在明顯變化。

為了進(jìn)一步分析各個腦區(qū)能量變化的規(guī)律，本文將大腦分為前額區(qū)(F)、左顳區(qū)(LT)、中心區(qū)(CP)、枕區(qū)(O)、右顳區(qū)(RT)五個腦區(qū)。接下來在各個腦區(qū)上，對癲癇發(fā)作前、中、后(R1，R2，R3)三個階段的EEG的功率譜密度(PSD)值進(jìn)行單因素方差分析，結(jié)果如圖7所示。當(dāng)p<0.01(圖中用“*”注明)，表明狀態(tài)之間存在顯著性差異。圖7中，在delta頻帶下(圖7(a))，腦區(qū)能量的變化與地形圖分析的結(jié)果一致，發(fā)作中各個腦區(qū)的能量值均高于發(fā)作前。與發(fā)作后相比，只有前額區(qū)、右顳區(qū)的能量值相近，并且發(fā)作期間各個腦區(qū)之間能量的變化沒有發(fā)作前、中顯著；發(fā)作中與發(fā)作前、后的能量只在右顳區(qū)沒有顯著性區(qū)別。delta頻帶下，在發(fā)作中各腦區(qū)能量的異常增強(qiáng)，而theta頻帶下(圖7(b))，除中心區(qū)以外，各腦區(qū)能量均低于發(fā)作前和發(fā)作后。從能量分布的區(qū)域性角度看，發(fā)作前、后不同腦區(qū)的能量存在較明顯的區(qū)別，而發(fā)作中各腦區(qū)能量值趨于一致；除中心區(qū)以外，其余各腦區(qū)發(fā)作中與前、后能量存在顯著差異。而alpha和beta頻帶下(圖7(c)、(d))的能量值較低，各腦區(qū)能量的變化無明顯規(guī)律。

圖7 4個頻帶下癲癇發(fā)作前、中、后各腦區(qū)功率譜密度的顯著性分析

接下來，對發(fā)作前、中、后能量存在顯著變化的theta頻帶下的第六通道EEG信號C4采用有限穿越可視圖(LPVG)進(jìn)行分析。序列長度選取為2 000個數(shù)據(jù)點(diǎn)。圖8為構(gòu)造的復(fù)雜網(wǎng)絡(luò)的鄰接矩陣，構(gòu)造的復(fù)雜網(wǎng)絡(luò)能較好地還原腦電序列的時域特性，并且對噪聲具有極強(qiáng)的適應(yīng)能力。由圖8所示可知，采用LPVG算法提取的發(fā)作前、中、后三種狀態(tài)下的EEG序列的鄰接矩陣中看出：theta頻帶下，矩陣中社團(tuán)連接主要集中在對角線上，其他區(qū)域的社團(tuán)連接面積都較小。發(fā)作前、后(圖8(a)、(c))，對角線上社團(tuán)連接大小隨機(jī)分布，而發(fā)作中(圖8(b))，社團(tuán)連接面積有所增大，但分布呈現(xiàn)一定的規(guī)律性，表明發(fā)作中網(wǎng)絡(luò)的復(fù)雜性有所降低。綜上，發(fā)作發(fā)作前、后，網(wǎng)絡(luò)模塊分布無規(guī)則，各模塊大小有明顯區(qū)分；而發(fā)作中，網(wǎng)絡(luò)中各模塊大小趨于一致，并且小規(guī)模的網(wǎng)絡(luò)各模塊明顯減少，各模塊的獨(dú)立性有所提高，網(wǎng)絡(luò)的復(fù)雜度降低，相對趨向于規(guī)則網(wǎng)絡(luò)。

圖8 有限穿越可視圖網(wǎng)絡(luò)矩陣

圖9為theta頻帶下網(wǎng)絡(luò)的節(jié)點(diǎn)度和度分布?？梢园l(fā)現(xiàn)，發(fā)作中(圖9(b))，網(wǎng)絡(luò)的度值明顯高于發(fā)作前、后(圖9(a)、(c))，且度分布的峰值增大且峰值右移，概率最大的度值接近40，而發(fā)作前、后概率最大的度值在20左右。而發(fā)作后與發(fā)作前相比，節(jié)點(diǎn)度分布要分散一些，但總體變化不明顯。

圖10為采用LPVG對癲癇發(fā)作前、中、后的EEG序列在delta、theta、alpha和beta四個頻帶下分別構(gòu)建，并提取網(wǎng)絡(luò)特征參數(shù)：平均路徑長度和圖指數(shù)復(fù)雜度(GIC)。由圖10(a)可知，在theta和alpha頻帶下，前、中、后三種狀態(tài)下區(qū)分不明顯；在delta頻帶下，發(fā)作中網(wǎng)絡(luò)的平均最短路徑高于發(fā)作前、后，表明網(wǎng)絡(luò)信息傳遞減慢；在beta頻帶下，發(fā)作中網(wǎng)絡(luò)的平均路徑長度變小，考慮到高頻帶包含的信息量較少，因此，對分析結(jié)果沒有明顯影響。圖10(b)可知，網(wǎng)絡(luò)的復(fù)雜度隨著頻帶增高而降低，在delta和theta頻帶下，發(fā)作中網(wǎng)絡(luò)的復(fù)雜度要低于前、后，其中在delta頻帶的衰減最為顯著，表明癲癇發(fā)作期間，雖然delta頻帶腦電能量急劇增強(qiáng)，但腦網(wǎng)絡(luò)的復(fù)雜度降低，與EEG信號呈周期性變化相對應(yīng)；而在alpha和beta頻帶下，復(fù)雜度改變不明顯，對應(yīng)高頻帶下網(wǎng)絡(luò)信息量較低。

圖9 theta頻帶下有限穿越可視圖網(wǎng)絡(luò)的節(jié)點(diǎn)度以及度值概率分布

圖10 4個頻帶下癲癇發(fā)作前、中、后的網(wǎng)絡(luò)特征

5 結(jié)語

為探究癲癇發(fā)作對腦電信號的影響，本文采用功率譜密度和有限穿越可視圖算法算法分析癲癇發(fā)作前、中、后三種狀態(tài)下的單通道腦電信號能量分布以及腦電信號的網(wǎng)絡(luò)特征。研究發(fā)現(xiàn)，癲癇發(fā)作時，腦電信號能量分布發(fā)生明顯改變，能量向低頻帶轉(zhuǎn)移，表現(xiàn)為delta頻帶能量顯著增加，而theta頻帶能量明顯降低。對腦區(qū)能量分布的進(jìn)一步研究發(fā)現(xiàn)，在delta頻帶下，發(fā)作前與發(fā)作后各腦區(qū)能量均有顯著性差異；在theta頻帶下，除中心區(qū)外，各腦區(qū)均有顯著性差異，表明癲癇發(fā)作對大部分腦區(qū)在低頻帶的腦電信號影響較大。對第六通道C4在theta頻帶下的發(fā)作前、中、后三種狀態(tài)下對應(yīng)的EEG信號進(jìn)行有限穿越可視圖(LPVG)分析，發(fā)現(xiàn)三種狀態(tài)下，LPVG構(gòu)建的復(fù)雜網(wǎng)絡(luò)具有不同形態(tài)的拓?fù)浣Y(jié)構(gòu)，發(fā)作中時更接近于規(guī)則網(wǎng)絡(luò)。與發(fā)作前、后相比，在高頻帶(alpha和beta)下，發(fā)作中網(wǎng)絡(luò)的平均最短路徑減小，而在低頻帶(delta和theta)顯著增大，且網(wǎng)絡(luò)的復(fù)雜度明顯降低。本文研究結(jié)果表明，功率譜密度和有限穿越可視圖算法能夠有效表征癲癇前、中、后三種狀態(tài)下腦電信號能量分布和單通道信號可視化后的網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)的異常，為癲癇的病理研究和臨床診斷提供幫助。

References)

[1] HANSEN B H, ALFSTAD K A, Van ROY B, et al.Sleep problems in children and adolescents with epilepsy: Associations with psychiatric comorbidity [J].Epilepsy & Behavior, 2016, 62: 14-19.

[2] PUKA K, RUBINGER L, CHAN C, et al.Predictors of intellectual functioning after epilepsy surgery in childhood: the role of socioeconomic status [J].Epilepsy & Behavior, 2016, 62: 35-39.

[3] NAVARRO-SUNE X, HUDSON A L, DE VICO FALLANI F, et al.Riemannian geometry applied to detection of respiratory states from EEG signals: the basis for a brain-ventilator interface [EB/OL].[2016-04-20].https://arxiv.org/pdf/1601.03022.

[4] SHEN C, LIU S, ZHOU W, et al.A physiology-based seizure detection system for multichannel EEG [J].PLoS One, 2013, 8(6): e65862.

[5] CZIGLER B, CSIKóS D, HIDASI Z, et al.Quantitative EEG in early Alzheimer’s disease patients-power spectrum and complexity features [J].International Journal of Psychophysiology, 2008, 68(1): 75-80.

[6] KLOVATCH-PODLIPSKY I, GAZIT T, FAHOUM F, et al.Dual array EEG-fMRI: an approach for motion artifact suppression in EEG recorded simultaneously with fMRI [EB/OL].[2016-03-20].http://www.sciencedirect.com/science/article/pii/S1053811916303214.

[7] SCHLEIGER E, WONG A, READ S, et al.Improved cerebral pathophysiology immediately following thrombectomy in acute ischaemic stroke: Monitoring via quantitative EEG [J].Clinical Neurophysiology : Official Journal of the International Federation of Clinical Neurophysiology, 2016, 127(8): 2832-2833.

[8] 李昌強(qiáng),黃力宇,鞠烽熾,等.癲癇腦電的雙譜特性研究[J].北京生物醫(yī)學(xué)工程,2004,23(1):9-13.(LI C Q, HUANG L Y, JU F Z, et al.Study of epileptic electroencephalogram using bispectrum analysis [J].Beijing Biomedical Engineering, 2004,23(1): 9-13.)

[9] 劉建平,張一聞.一種高分辨率時頻分析方法及在腦電信號處理中的應(yīng)用[J].武警工程學(xué)院學(xué)報,2002,18(6):19-22.(LIU J P, ZHANG Y W.A high resolution time frequency analysis method and its application in EEG signal processing [J].Journal of Engineering University of the Chinese People’s Armed Police Force, 2002, 18(6): 19-22.)

[10] ADELI H, ZHOU Z, DADMEHR N.Analysis of EEG records in an epileptic patient using wavelet transform [J].Journal of Neuroscience Methods, 2003, 123(1): 69-87.

[11] MPORAS I, TSIRKA V, ZACHARAKI E I, et al.Seizure detection using EEG and ECG signals for computer-based monitoring, analysis and management of epileptic patients [J].Expert Systems with Applications, 2015, 42(6): 3227-3233.

[12] 杜守洪,陳子怡,田翔華,等.近似熵與SVM在自動分類癲癇腦電信號中的應(yīng)用[J].生物醫(yī)學(xué)工程研究,2013,32(2):74-79.(DU S H, CHEN Z Y, TIAN X H, et al.Epilepsy EEG Approximate Entropy (ApEn) Support Vector Machine (SVM) real-time detection classification [J].Journal of Biomedical Engineering Research, 2013, 32(2):74-79.)

[13] 張濤,陳萬忠,李明陽.基于AdaBoost算法的癲癇腦電信號識別[J].物理學(xué)報,2015,64(12):128701.(ZHANG T, CHEN W Z, LI M Y.Recognition of epilepsy electroencephalography based on AdaBoost algorithm [J].Acta Physica Sinica, 2015, 64(12): 128701.)

[14] XIANG J, LI C, LI H, et al.The detection of epileptic seizure signals based on fuzzy entropy [J].Journal of Neuroscience Methods, 2015, 243:18-25.

[15] TAKAHASHI T, CHO R Y, MIZUNO T, et al.Antipsychotics reverse abnormal EEG complexity in drug-naive schizophrenia: a multiscale entropy analysis [J].Neuroimage, 2010, 51: 173-182.

[16] SONG Y, CROWCROFT J, ZHANG J.Automatic epileptic seizure detection in EEGs based on optimized 21 sample entropy and extreme learning machine [J].Journal of Neuroscience Methods, 2012, 210: 132-146.

[17] FELDHOFF J H, DONNER R V, DONGES J F, et al.Geometric detection of coupling directions by means of inter-system recurrence networks [J].Physics Letters A, 2013, 376(46): 3504-3513.

[18] DEZA J I, BARREIRO M, MASOLLER C.Inferring interdependencies in climate networks constructed at inter-annual, intra-season and longer time scales [J].European Physical Journal Special Topics, 2013, 222(2):511-523.

[19] FELDHOFF J H, DONNER R V, DONGES J F, et al.Geometric detection of coupling directions by means of inter-system recurrence networks [J].Physics Letters A, 2012, 376(46):3504-3513.

[20] ZHANG J, SMALL M.Complex network from pseudoperiodic time series: topology versus dynamics [J].Physical Review Letters, 2006, 96(23): 238701-238701.

[21] GAO Z, JIN N.Flow-pattern identification and nonlinear dynamics of gas-liquid two-phase flow in complex networks [J].Physical Review E: Statistical Nonlinear & Soft Matter Physics, 2009, 79(6): 1019-1027.

[22] WANG R, WANG J, YU H, et al.Decreased coherence and functional connectivity of electroencephalograph in Alzheimer’s disease [J].Chaos, 2014, 24(3): 033136-033136.

[23] LACASA L, LUQUE B, BALLESTEROS F.From time series to complex networks: the visibility graph [J].Proceedings of the National Academy of Sciences of the United States of America, 2008, 105(13): 4972-4975.

[24] MEHRAN A, HOJJAT A.Visibility graph similarity: a new measure of generalized synchronization in coupled dynamic systems [J].Physica D: Nonlinear Phenomena, 2012, 241(4): 326-332.

[25] MEHRAN A, HOJJAT A, AMIR A.Improved visibility graph fractality with application for the diagnosis of autism spectrum disorder [J].Physica A: Statistical Mechanics & Its Applications, 2012, 391(20): 4720-4726.

[26] JIANG S, BIAN C, NING X, et al.Visibility graph analysis on heartbeat dynamics of meditation training [J].Applied Physics Letters, 2013, 102(25): 253702-1-253702-3.

[27] TELESCA L, LOVALLO M, PIERINI J O.Visibility graph ap-proach to the analysis of ocean tidal records [J].Chaos Solitons & Fractals, 2012, 45(9/10): 1086-1091.

[28] CHARAKOPOULOS A K, KARAKASIDIS T E, PAPANICOLAOU P N, et al.The application of complex network time series analysis in turbulent heated jets [J].Chaos, 2014, 24(2): 024408-024408.

[29] WANG N, LI D, WANG Q.Visibility graph analysis on quarterly macroeconomic series of China based on complex network theory [J].Physica A: Statistical Mechanics & Its Applications, 2012, 391(24): 6543-6555.

[30] LUQUE B, LACASA L, BALLESTEROS F, et al.Horizontal visibility graphs: exact results for random time series [J].Physical Review E: Statistical Nonlinear & Soft Matter Physics, 2009, 80(2): 593-598.

[31] 周婷婷,金寧德,高忠科,等.基于有限穿越可視圖的時間序列網(wǎng)絡(luò)模型[J].物理學(xué)報,2012,61(3):030506.(ZHOU T T, JIN N D, GAO Z K, et al.Limited penetrable visibility graph for establishing complex network from time series [J].Acta Physica Sinica, 2012, 61(3): 030506.)

[32] PEI X, WANG J, DENG B, et al.WLPVG approach to the analysis of EEG-based functional brain network under manual acupuncture [J].Cognitive Neurodynamics, 2014, 8(5): 417-428.

[33] AKIN M, KIYMIK M K.Application of periodogram and AR spectral analysis to EEG signals [J].Journal of Medical Systems, 2000, 24(4): 247-256.

[34] LI L, XIAO L, CHEN L.Differences of EEG between eyes-open and eyes-closed states based on autoregressive method [J].Journal of Electronic Science and Technology of China, 2009, 7(2): 175-179.

[35] STAM C J, REIJNEVELD J C.Graph theoretical analysis of complex networks in the brain [J].Nonlinear Biomedical Physics, 2007, 1: 3.

[36] KIM J, WILHELM T.What is a complex graph? [J].Physica A: Statistical Mechanics and its Applications, 2008, 387(11): 2637-2652.

[37] BOCCALETTI S, LATORA V, MORENO Y, et al.Complex networks: structure and dynamics [J].Complex Systems & Complexity Science, 2007, 424(4/5): 175-308.

This work is partially supported by the National Natural Science Foundation of China (61302002, 61601331), the Natural Science Foundation of Tianjin City (14JCYBJC15400, 14JCQNJC01200), Tangshan Technology Research and Development Program (14130223B), Advance Research Program of Tianjin University of Technology and Education (KYQD1611).

WANG Ruofan, born in 1986, Ph.D., lecturer.Her research interests include data analysis and algorithm mining, analysis of complex nonlinear system.

LIU Jing, born in 1973, Ph.D., chief physician.Her research interests include neurology.

WANG Jiang, born in 1964, Ph.D., professor.His research interests include neural computing.

YU Haitao, born in 1985, Ph.D., lecturer.His research interests include nervous system modeling and kinetic analysis.

CAO Yibin, born in 1963, M.S., chief physician.His research interests include neurology.

Analysis algorithm of electroencephalogram signals for epilepsy diagnosis based on power spectral density and limited penetrable visibility graph

WANG Ruofan1, LIU Jing2, WANG Jiang3, YU Haitao3, CAO Yibin2*

(1.CollegeofInformationTechnologyEngineering,TianjinUniversityofTechnologyandEducation,Tianjin300222,China;2.DepartmentofNeurology,TangshanGongrenHospital,TangshanHebei064300,China;3.CollegeofElectrical&AutomationEngineering,TianjinUniversity,Tianjin300072,China)

Focused on poor robustness to noise of the Visibility Graph (VG) algorithm, an improved Limited Penetrable Visibility Graph (LPVG) algorithm was proposed.LPVG algorithm could map time series into networks by connecting the points of time series which satisfy the certain conditions based on the visibility criterion and the limited penetrable distance.Firstly, the performance of LPVG algorithm was analyzed.Secondly, LPVG algorithm was combined with Power Spectrum Density (PSD) to apply to the automatic identification of epileptic ElectroEncephaloGram (EEG) before, during and after the seizure.Finally, the characteristic parameters of the LPVG network in the three states were extracted to study the influence of epilepsy seizures on the network topology.The simulation results show that compared with VG and Horizontal Visibility Graph (HVG), although LPVG had a high time complexity, it had strong robustness to noise in the signal: when mapping the typical periodic, random, fractal and chaos time series into networks by LPVG, it was found that as the noise intensity increased, the fluctuation rates of clustering coefficient by LPVG network were always the lowest, respectively 6.73%, 0.05%, 0.99% and 3.20%.By the PSD and LPVG analysis, it was found that epilepsy seizure had great influence on the brain energy.PSD was obviously enhanced in the delta frequency band, and significantly reduced in the theta frequency band; the topological structure of the LPVG network changed during the seizure, characterized by the independent enhanced network module, increased average path length and decreased graph index complexity.The PSD and LPVG applied in this paper could be taken as an effective measure to characterize the abnormality of the energy distribution and topological structure of single EEG signal channel, which would provide help for the pathological study and clinical diagnosis of epilepsy.

ElectroEncephaloGram (EEG); epilepsy; Power Spectral Density (PSD); Limited Penetrable Visibility Graph (LPVG); complex network

2016-08-10;

2016-09-28。 基金項(xiàng)目:國家自然科學(xué)基金資助項(xiàng)目(61302002， 61601331)；天津市自然科學(xué)基金資助項(xiàng)目(14JCYBJC15400， 14JCQNJC01200)；唐山市科技支撐項(xiàng)目(14130223B)；天津職業(yè)技術(shù)師范大學(xué)預(yù)研項(xiàng)目(KYQD1611)。

王若凡(1986—)，女，河北滄州人，講師，博士，主要研究方向：數(shù)據(jù)分析與算法挖掘、復(fù)雜非線性系統(tǒng)分析； 劉靜(1973—)，女，河北唐山人，主任醫(yī)師，博士，主要研究方向：神經(jīng)病學(xué)； 王江(1964—)，男，河北唐山人，教授，博士，主要研究方向：神經(jīng)計算； 于海濤(1985—)，男，河北唐山人，講師，博士，主要研究方向：神經(jīng)系統(tǒng)建模與動力學(xué)分析； 曹亦賓(1963—)，男，河北唐山人，主任醫(yī)師，碩士，主要研究方向：神經(jīng)病學(xué)。

1001-9081(2017)01-0175-08

10.11772/j.issn.1001-9081.2017.01.0175

TP391.4

A