朱寧寧,李皓,鄧小喬,于明,李效龍
1.江蘇科技大學(xué)電子信息學(xué)院,江蘇鎮(zhèn)江212000;2.江蘇大學(xué)附屬醫(yī)院神經(jīng)內(nèi)科,江蘇鎮(zhèn)江212000
癲癇是多種病因引起的慢性腦功能障礙綜合征,主要是由大腦神經(jīng)細(xì)胞群反復(fù)超同步放電而引起的,具有突發(fā)性和短暫性腦功能障礙等特點(diǎn),會(huì)引起患者認(rèn)知功能損害和情緒障礙。癲癇患者腦電信號(hào)中會(huì)出現(xiàn)棘波、尖波、棘慢復(fù)合波、尖慢復(fù)合波等異常波形。由于棘波幅度大、電位短、垂直上升和下降等瞬態(tài)特征,使得癲癇腦電信號(hào)中棘波的檢測(cè)具有臨床診斷意義。常見(jiàn)的棘波波形多為負(fù)相棘波,有時(shí)也為正相棘波,還有雙相、三相棘波等,具體波形如圖1所示[1]。
目前常用的癲癇發(fā)作檢測(cè)算法有形態(tài)成分分析[2]、獨(dú)立分量分析[3]、小波變換[4]、希爾伯特黃變換[5]、頻率小波變換和支持向量機(jī)[6]、遞歸量化分析與支持向量機(jī)[7]、非矩陣分解[8]、奇異值分解[9]、深度卷積網(wǎng)絡(luò)[10]、改進(jìn)遺傳算法[11]、稀疏表示[12-16]等。
圖1 不同形態(tài)結(jié)構(gòu)的棘波Fig.1 Spikes of different morphological structures
癲癇發(fā)作檢測(cè)要求腦電信號(hào)在特定的時(shí)間段或頻率段有盡可能高的頻率或時(shí)間分辨率,同時(shí)要求較高的棘波檢出率以提高檢測(cè)準(zhǔn)確率。本文結(jié)合小波包變換可以多層次劃分信號(hào)頻帶,對(duì)癲癇腦電信號(hào)進(jìn)行3層小波包分解得到每層的腦電波頻率帶,然后重構(gòu)腦電頻率范圍內(nèi)關(guān)鍵節(jié)點(diǎn)頻率的癲癇腦電信號(hào),最后再結(jié)合棘波的幅度特征提取出癲癇患者在不同時(shí)期的棘波。
小波包變換不僅可以很好地表征以低頻信息為主要成分的癲癇腦電信號(hào),還能夠?qū)π盘?hào)高頻部分進(jìn)行更精細(xì)的分解,并且可以依據(jù)信號(hào)特點(diǎn)自適應(yīng)的選擇對(duì)應(yīng)的頻帶,使得信號(hào)的頻帶與頻譜相匹配,從而提高信號(hào)的時(shí)頻分辨率,因此是一種適用于癲癇腦電棘波檢測(cè)的信號(hào)分析方法。采用小波包變換的棘波檢測(cè)算法流程如圖2所示。
1.1.1 癲癇腦電信號(hào)頻帶劃分實(shí)驗(yàn)數(shù)據(jù)的采樣頻率為173.61 Hz,信號(hào)長(zhǎng)度為23.6 s。根據(jù)奈奎斯特采樣定理可知,信號(hào)的采樣頻率范圍為0~86.80 Hz。對(duì)采集的癲癇腦電信號(hào)進(jìn)行三層小波包分解得到腦電信號(hào)的頻帶劃分,具體如表1所示。
1.1.2 重構(gòu)特定節(jié)點(diǎn)頻率的腦電信號(hào)從表1 可以看出,第三層小波包分解把采集的癲癇信號(hào)進(jìn)行了更細(xì)致的分解。由于大部分腦電信號(hào)的頻率范圍在0~30 Hz,因?yàn)楸? 中節(jié)點(diǎn)S(3,0)、S(3,1)、S(3,2)包含了所有腦電波的頻率,也包括要提取的棘波頻率,因此只需重構(gòu)節(jié)點(diǎn)S(3,0)、S(3,1)、S(3,2)即能夠最大程度上得到所要提取的癲癇患者棘波信號(hào)。
三層小波包分解如圖3所示,小波包變換中的每一層信號(hào)都被分解成低頻和高頻部分[17]。例如信號(hào)S經(jīng)過(guò)三層小波包分解之后可以表示為:
圖2 棘波檢測(cè)流程圖Fig.2 Flowchart of epileptic spike detection
表1 基于小波包分解的腦電信號(hào)頻帶劃分Tab.1 Electroencephalogram(EEG)signal frequency division based on wavelet packet decomposition
可見(jiàn),信號(hào)S經(jīng)過(guò)三層小波包分解后的子頻帶包含了原信號(hào)所有的細(xì)節(jié)頻率信息。
圖3 三層小波包分解樹(shù)結(jié)構(gòu)圖Fig.3 Diagram of 3-layer wavelet packet decomposition tree
1.1.3 棘波檢測(cè)閾值選取棘波具有的明顯的物理特征[18],其中時(shí)限為1/50~1/14 s,放入腦電圖測(cè)量尺14 Hz的刻度內(nèi)則為20~70 ms,幅度多為100~200 μV。本文采用小波包變換對(duì)節(jié)點(diǎn)S(3,0)、S(3,1)、S(3,2)的腦電信號(hào)進(jìn)行重構(gòu)后,選擇棘波信號(hào)的幅度作為閾值進(jìn)行癲癇腦電棘波的提取。
對(duì)腦電信號(hào)f(t) ∈L2(R),設(shè)它的小波分解為:
fj(t) ∈Wj的小波級(jí)數(shù)為:
根據(jù)文獻(xiàn)[19]給出的定理,對(duì)于每個(gè)j=1,2,…有:
進(jìn)而,對(duì)于每個(gè)m=0,1,…,2k -1,k=1,2,…,j和j=2,…函數(shù)族,有:
對(duì)于正交小波包,根據(jù)文獻(xiàn)[19]給出的定理,?k(1 ≤k≤j)和m(0 ≤m≤2k- 1),fj(t)∈Wj可以進(jìn)一步分解成小波包分量的正交和,即:
令k=0,則由式(6)可得f jn(t)∈U jn,式(5)可改寫為:
正交小波基定義:令ωn(t)滿足下列尺度方程:
其中,n=0,1,2,…,系數(shù)序列{hn} 和{gn} 仍滿足下列關(guān)系:
其中,h、g為濾波器系數(shù),d為小波包分解系數(shù),l、k為分解層數(shù),j、n為小波包節(jié)點(diǎn)號(hào)。
實(shí)驗(yàn)數(shù)據(jù)來(lái)自德國(guó)Born 癲癇腦電研究室。該數(shù)據(jù)庫(kù)中共包含標(biāo)號(hào)從A~E 的5 個(gè)數(shù)據(jù)集(表2),其中A、B 數(shù)據(jù)集是來(lái)自健康志愿者的腦電信號(hào),C、D、E的數(shù)據(jù)集是癲癇患者手術(shù)前診斷并記錄的顱內(nèi)腦電信號(hào)。C、D 是癲癇發(fā)作間期采集的數(shù)據(jù),E 是癲癇發(fā)作期采集的信號(hào)。每個(gè)數(shù)據(jù)集均包括100個(gè)采用10-20國(guó)際標(biāo)準(zhǔn)導(dǎo)聯(lián)技術(shù)提取的單導(dǎo)聯(lián)腦電信號(hào)。每個(gè)腦電信號(hào)的長(zhǎng)度為23.6 s,采樣頻率為173.6 Hz。數(shù)據(jù)經(jīng)過(guò)了一定程度的預(yù)處理,去除了明顯的噪聲和偽跡。仿真采用的數(shù)據(jù)包括患者癲癇發(fā)作期的數(shù)據(jù)、健康期的數(shù)據(jù)和癲癇發(fā)作間期的數(shù)據(jù)。實(shí)驗(yàn)的仿真平臺(tái)基于Matlab2016b。
表2 德國(guó)Born癲癇腦電研究室數(shù)據(jù)Tab.2 Epileptic EEG data from Born epilepsy laboratory,Germany
對(duì)癲癇患者a的腦電信號(hào)進(jìn)行頻譜分析(圖4)和三層小波包分解(圖5),然后采用本文提出的算法對(duì)癲癇患者a 在健康期、癲癇發(fā)作間期的腦電信號(hào)進(jìn)行棘波提取,分別如圖6 和圖7 所示。因?yàn)榘d癇發(fā)作期棘波出現(xiàn)的比較多,為了更準(zhǔn)確的提取在腦電頻率段出現(xiàn)的棘波,根據(jù)表1的腦電頻帶劃分將癲癇發(fā)作期棘波的提取劃分為節(jié)點(diǎn)頻率(3, 0)的提取和節(jié)點(diǎn)頻率(3, 1)的提取。為了檢測(cè)該算法是否具有特異性(該算法是否對(duì)其他癲癇患者同樣有效),分別將癲癇患者a 和癲癇患者b 癲癇發(fā)作期的數(shù)據(jù)進(jìn)行對(duì)比,以驗(yàn)證該算法不存在特異性,分別如圖8、圖9、圖10、圖11所示。
最后根據(jù)仿真提取棘波的檢測(cè)結(jié)果總結(jié)了5 個(gè)癲癇患者在不同時(shí)期的棘波檢測(cè)數(shù)、誤檢數(shù)、漏檢數(shù)以及相對(duì)應(yīng)的各個(gè)時(shí)期的誤檢率和漏檢率(表3)。
仿真結(jié)果分析如下所示。
圖4 癲癇患者a健康期的腦電信號(hào)Fig.4 EEG signals of patient a in healthy period
圖5 小波包分解樹(shù)結(jié)構(gòu)及節(jié)點(diǎn)系數(shù)Fig.5 Wavelet packet decomposition tree and node coefficients
(1)由于仿真的數(shù)據(jù)已經(jīng)過(guò)預(yù)處理,去除了噪聲和偽跡,因此利用小波包變換和棘波的頻率、幅度等物理特征相結(jié)合的算法不僅能提取出同一個(gè)病人在健康期、癲癇發(fā)作間期、癲癇發(fā)作期等不同時(shí)期的棘波;而且可以提取出不同患者在癲癇發(fā)作期的棘波,證明該算法適用于不同癲癇患者的棘波檢測(cè),不存在特異性。
(2)癲癇患者腦電中的棘波主要以爆發(fā)式節(jié)律出現(xiàn)。由健康期、癲癇發(fā)作間期、癲癇發(fā)作期相應(yīng)節(jié)點(diǎn)檢測(cè)出的棘波可以發(fā)現(xiàn)棘波主要出現(xiàn)在0~22 Hz頻率段。所提取出的棘波波形包含參考文獻(xiàn)[14]所指出的14 Hz 或6 Hz 的正相棘波,當(dāng)然這里的6 Hz通常指5~7 Hz 也包含13~17 Hz,而且各自以獨(dú)個(gè)連續(xù)爆發(fā)的形式出現(xiàn)。
圖6 癲癇患者a健康期腦電信號(hào)處理Fig.6 EEG signal processing of patient a in healthy period
圖7 癲癇患者a發(fā)作間期腦電信號(hào)處理Fig.7 EEG signal processing of patient a in intermittent period
圖8 癲癇患者a發(fā)作期節(jié)點(diǎn)(3,0)的腦電信號(hào)處理Fig.8 EEG signal processing of patient a in node(3,0)during epilepsy attack
圖9 癲癇患者b發(fā)作期節(jié)點(diǎn)(3,0)的腦電信號(hào)處理Fig.9 EEG signal processing of patient b in node(3,0)during epilepsy attack
圖10 癲癇患者a發(fā)作期節(jié)點(diǎn)(3,1)的腦電信號(hào)處理Fig.10 EEG signal processing of patient a in node(3,1)during epilepsy attack
圖11 癲癇患者b發(fā)作期節(jié)點(diǎn)(3,1)的腦電信號(hào)處理Fig.11 EEG signal processing of patient b in node(3,1)during epilepsy attack
表3 5個(gè)癲癇患者在不同時(shí)期的棘波檢測(cè)結(jié)果Tab.3 Results of spike detection in 5 epilepsy patients at different time
(3)基于小波包變換的棘波提取算法在癲癇患者健康期和癲癇發(fā)作間期的棘波提取準(zhǔn)確率都在85%以上。癲癇患者發(fā)作期的平均誤檢率為14.56%(誤檢率=誤檢數(shù)/專家檢測(cè)數(shù),平均誤檢率為癲癇發(fā)作期誤檢率的平均值);平均漏檢率為10.57%(漏檢率=漏檢數(shù)/專家檢測(cè)數(shù),平均漏檢率為癲癇發(fā)作期漏檢率的平均值)。整體癲癇患者的棘波誤檢率為12.02%(總誤檢數(shù)/總專家檢測(cè)數(shù));整體漏檢率為11.70%(總漏檢數(shù)/總專家檢測(cè)數(shù))。癲癇病人發(fā)作期棘波的誤檢率、漏檢率差別比較大是因?yàn)橛械募〝?shù)據(jù)在腦電圖上不是很明顯,在人為檢測(cè)時(shí)會(huì)造成一定程度的誤檢、漏檢,而且選取的棘波物理特征不全面,例如沒(méi)有考慮棘波的脈寬這一特性也使得該算法會(huì)造成一定程度的漏檢和誤檢。
本文提出的基于小波包變換結(jié)合棘波的頻率、幅度等物理特性提取棘波的算法,能夠提取出癲癇病人發(fā)作期腦電信號(hào)里大部分棘波。仿真結(jié)果得出該算法具有很低的漏檢率和誤檢率(分別為12.02%和11.70%)。分析得出的結(jié)論與參考文獻(xiàn)[14]的相關(guān)結(jié)論吻合。該算法可以為醫(yī)生判斷癲癇是否發(fā)作,進(jìn)而確定癲癇的病灶提供了參考。