腦電信號的最優(yōu)分數(shù)階傅里葉變換

張相芬, 陀佳萍, 董柳吟, 袁非牛, 羅 陽

(上海師范大學 信息與機電工程學院, 上海200234)

腦電信號是頭皮或大腦皮層腦神經(jīng)細胞生理活動的反應, 它包含了許多病理信息, 可以為某些腦部疾病提供診斷依據(jù)[1]. 自19世紀20年代檢測到腦電信號以來, 人們對其開展了大量的研究工作, 腦電信號的處理和分析至今仍然是一項十分困難但又非常重要的研究課題[2].

腦電信號的處理主要分為2個步驟[3]:首先是去除眼電偽跡、去直流和干擾,即將非腦電信號成分去除的預處理;其次是對腦電信號進行特征提取、模式識別等分析,以進行診斷分類或者識別.目前,腦電信號的預處理方法主要包括表面拉普拉斯參考法、主成分分析法及獨立成分分析法(independent component analysis,ICA)等.表面拉普拉斯參考法可以在一定程度上去除冗余信號,但未考慮不同位置頭皮傳輸特性的差異[4];主成分分析法可用于移除多通道腦電信號中的眼電干擾成分,此方法的缺點是去線性相關的約束不夠[5];ICA主要是去除信號的獨立成分干擾[6].

常用的腦電信號處理方法包括非線性動力學分析、人工神經(jīng)網(wǎng)絡分析、高階譜分析以及時頻分析等.隨著非線性學科的發(fā)展,非線性動力學成為分析和處理腦電信號的有效工具.但是,非線性動力學分析方法對初值極其敏感,因此采用非線性動力學分析處理腦電信號存在困難[7].人工神經(jīng)網(wǎng)絡能高速尋找最優(yōu)化解,但其泛化能力比較弱,較難找到通用模型,且容易陷入局部最優(yōu)解[8].高階譜分析在處理非線性腦電信號時能夠有效抑制高斯噪聲,但在真值附近起伏很大,受噪聲干擾影響嚴重[9].

腦電信號較為復雜,從定量分析的角度來看,采用傳統(tǒng)的時域或者頻域分析方法往往性能不佳,因此,研究者想到采用時域與頻域結合的方法,即時頻分析方法.分數(shù)階傅里葉變換是一種典型的時頻分析方法,在腦電信號分析中受到了關注[10].腦電信號在最優(yōu)分數(shù)階傅里葉變換域有很好的能量聚集性[11].本文基于二維峰值搜索算法提出基于最優(yōu)分數(shù)階傅里葉變換對腦電信號進行分析處理,以為腦電信號特征提取及分類提供一種新的途徑.

1 腦電信號及其預處理

1.1 腦電信號

腦電信號非常微弱,正常情況下,腦電信號的幅值在5～100 μV之間,其頻率范圍為0.5～35.0 Hz[12].按照頻帶劃分,可將腦電信號分為4類[13].

1)δ波,極低的頻帶:0.5～3.0 Hz,幅值:10～20 μV,出現(xiàn)在睡眠中、人體缺氧或者大腦病變時.

2)θ波,低頻帶:4～7 Hz,幅值:20～40 μV,正常成年人在清醒狀態(tài)下不會表現(xiàn)此波,是一種與幻想、睡眠或者回憶狀態(tài)有關的波段信號.

3)α波,中頻帶:8～13 Hz,幅值:10～100 μV,是正常人腦電信號的基本節(jié)律,普通人的α波節(jié)律調節(jié)較好,幅度和頻率變化都不明顯.

4)β波,高頻帶:14～26 Hz,幅值:5～30 μV,當個體處于情緒激動、緊張或者亢奮的狀態(tài)時,β波出現(xiàn)較多.

1.2 腦電信號預處理

腦電信號的預處理主要包括以下步驟[14].

1) 信號濾波.信號濾波即對原始的腦電信號進行低通及高通濾波[15].本文采用的腦電信號為MATLAB R2016a的EEGLAB工具箱自帶的數(shù)據(jù)[16].腦電信號的頻率范圍一般在0.5～35.0 Hz之間,因此濾波時先對信號進行截止頻率為0.5 Hz的高通濾波,再對其進行截止頻率為35 Hz的低通濾波.

2) 獨立成分分析.獨立成分分析是從多元(多維)統(tǒng)計數(shù)據(jù)中尋找潛在因子或成分的一種方法.將獨立成分分析應用在腦電信號分析中,可以較好地識別并去掉眼動和其他噪聲[17].

3) 眼電偽跡去除.利用Adjust工具進行眼電偽跡半自動去除[18].去除眼電偽跡后的信號波形更加平穩(wěn).

4) 疊加平均.即把腦電信號的每一段32個通道的數(shù)據(jù)進行疊加平均.

2 分數(shù)階傅里葉變換及二維峰值搜索算法

分數(shù)階傅里葉變換是傳統(tǒng)傅里葉變換的一種廣義形式.分數(shù)階傅里葉變換反映了信號在時域、頻域的信息,適用于非平穩(wěn)信號的處理,它的快速離散算法較為成熟,計算速度快[19].

分數(shù)階傅里葉變換可以有若干種不同的定義,目前常用的是Ozakats從積分變換角度給出的定義[20],即在時間域的函數(shù)x(t)的p階分數(shù)階傅里葉變換是一個線性積分運算,

(1)

其中,Rp是分數(shù)階傅里葉變換算子,核函數(shù)kp(t,v)的表達式為

(2)

式(1)可看作信號x(t)的p階分數(shù)階傅里葉變換,也可看作信號x(t)在α角度下的分數(shù)階傅里葉變換.因此,分數(shù)階傅里葉變換算子可記為Rα,其核函數(shù)也可記為kα(t,v).

本文采用的數(shù)值計算方法是Ozaktas采樣型離散化算法[21],即首先對原始信號進行時域解調,然后利用香農公式對解調后的信號進行插值計算,最后,對分數(shù)階域變量做離散化處理.

二維峰值搜索算法可以實現(xiàn)線性調頻信號的檢測與參數(shù)估計[22].基于二維峰值搜索算法尋找腦電信號的分數(shù)階傅里葉變換的最優(yōu)分數(shù)階次可概括為:以變換階次p作為變量進行掃描,對腦電信號進行連續(xù)的分數(shù)階傅里葉變換,形成信號在(p,v)平面上的二維分布.在(p,v)平面進行閾值做峰值點的二維搜索,從而實現(xiàn)腦電信號的最優(yōu)分數(shù)階次參數(shù)估計[23].

設信號的時域表達式為

其中:A為腦電信號的幅值;φ0為初始相位,這里設為0;f0為起始頻率;k0為腦電信號的調頻率,其值為信號帶寬與脈沖寬度的比.

信號的分數(shù)階傅里葉變換為

其中,

(5)

滿足k0+cot(pπ/2)=0的p為分數(shù)階變換的最優(yōu)變換階次,記為p0,并且有

3 結果與分析

我們從腦電信號(共80多個數(shù)據(jù)段)中挑選出第3(epochs 3)、第7(epochs 7)、第20(epochs 20)、第33(epochs 33)和第52(epochs 52)這5個數(shù)據(jù)段.疊加平均后的各個數(shù)據(jù)段信號時域波形分別如圖1所示.將疊加平均后的數(shù)據(jù)即腦電信號的幅值代入式(3),得到腦電信號的時域表達式.

圖1 各個數(shù)據(jù)段疊加平均后的時域波形Fig.1 Time domain waveforms after each segment is superimposed on average

在參數(shù)(p,v)平面上,以0.001的步進值對疊加平均后的信號在0

圖2為上述5個波段信號的二維峰值搜索及其一維峰值搜索的結果.由圖2可見,二維搜索能得出唯一的峰值,而一維搜索有可能出現(xiàn)多個峰值,較難得出最優(yōu)分數(shù)階.

圖2 各波段信號一維和二維峰值搜索對比圖Fig.2 Peak search contrast map of 1D and 2D signals in different bands

對選取的各個數(shù)據(jù)段的腦電信號進行了最優(yōu)分數(shù)階傅里葉變換與普通階次傅里葉變換,其結果如圖3所示.由圖3可以看出,與進行普通階次變換的腦電信號相比,進行最優(yōu)分數(shù)階變換的腦電信號的時頻域波形較為平穩(wěn),表明噪聲得到了抑制.

基于最優(yōu)分數(shù)階傅里葉變換處理后的腦電信號可以進行大腦的生理狀態(tài)判別.例如,圖3a所示為ephochs 3進行最優(yōu)分數(shù)階傅里葉變換的時頻域波形,顯然,0～0.8 s之間的波形為θ波,表明此時個體處于睡眠或者回憶狀態(tài);0.8～1.2 s之間的波形為α波,表明此時個體處于安靜或者閉眼狀態(tài);1.2～1.8 s之間的波形為β波,表明個體處于較激動亢奮的狀態(tài).而圖3b所示的普通階次變換下的腦電波形較難提取其節(jié)律波.

圖3 腦電信號最優(yōu)分數(shù)階變換與普通階次變換波形對比圖

由此可見,本文基于最優(yōu)分數(shù)階傅里葉變換的腦電信號處理對后續(xù)腦電信號的研究具有重要的意義.

4 結 語

腦電信號包含了大量的生理和病理信息,腦電信號可以為某些腦疾病的臨床診斷和治療提供依據(jù).本文的研究著眼于利用二維峰值搜索算法以及最優(yōu)分數(shù)階傅里葉變換對腦電信號進行處理和分析.

與一維搜索算法相比,二維搜索算法計算得到的最優(yōu)階次分數(shù)階傅里葉變換具有更好的能量聚集性.對比普通階次傅里葉變換的結果可知,基于二維峰值搜索算法進行腦電信號的最優(yōu)分數(shù)階傅里葉變換,能更好地去除信號的噪聲,讓信號波形更加平穩(wěn),為腦電信號的進一步研究打下了良好的基礎.

本文基于二維峰值搜索算法確定分數(shù)階傅里葉變換的最優(yōu)階次,并基于最優(yōu)階次的分數(shù)階傅里葉變換對腦電信號進行分析.

本研究采用的是EEGLAB自帶的腦電信號,信號成分較好,處理相對容易,因此該方法的適應性仍需進一步驗證;在使用二維峰值搜索算法尋找變換階次時,會存在計算速度與精度的矛盾.因此,最優(yōu)階次的搜索算法仍需進一步加以改進.