基于大鼠腦電信號時空特征精細回歸樹模型的步態(tài)解碼方法

寇建閣,李想,王娜,汪首坤,宮妍竹

（1.北京航空航天大學(xué)，北京 100191；2.北京理工大學(xué)，北京 100081；3.中國傳媒大學(xué)媒體融合與傳播國家重點實驗室，北京 100024）

1 引言

隨著腦-機技術(shù)的發(fā)展，獲取的腦電信號可靠性越來越高，足以作為控制信號來驅(qū)動外部設(shè)備［1,2］。但腦電信號的處理困難一直限制著腦-機技術(shù)的發(fā)展與工程實現(xiàn)，降低了腦電信號的實用性。不同的濾波方式、信號處理方法都會影響腦電信號的解碼質(zhì)量。

腦電信號一般有高通濾波和低通濾波的預(yù)處理方法。高通濾波器一般會得到高頻率峰電位（Spikes）信號；低通濾波后可以得到相對頻率相對較低的局部場電位信號（Local Potential Field,LFP)，此類腦電信號中蘊含著豐富的運動相關(guān)信息［3］。為了探究癲癇的成病因素，Gotman等人［4］利用截止頻率為1 Hz的2階低頻濾波器對腦電信號做了預(yù)處理，然后提取了預(yù)處理之后的腦電信號特征，對癲癇疾病患者的腦電信號判斷起到了很好的效果，有利于癲癇疾病的診斷。Hammer［5］、Takufumi等人［6］分別證明了腦電信號在解析運動中的價值，并通過實驗加以佐證。Jeong-Hun［7］等應(yīng)用低頻皮層腦電信號進行運動軌跡解碼，得出“低頻信號中蘊含了大量與運動有關(guān)的信息，包括運動開始時間、方位和速度等”的結(jié)論。

腦電信號濾波作為信號解碼的前提條件，取得了積極的成效，但腦電信號在實際分析應(yīng)用中仍會存在不同的問題。例如尖峰信號的分析中容易受到場地環(huán)境干擾，且長期記錄效果不佳，LFP信號分析中低頻動作、工頻干擾影響較大。因此，如何獲得一個穩(wěn)定、低噪的腦電信號是一重要關(guān)鍵問題。

獲得高質(zhì)量的預(yù)處理腦電信號之后，對該信號進行準確解碼是腦機接口的另一核心問題。目前，支持向量機（Support Vector Machines,SVM）、線性判別分析（Linear Discriminant Analysis,LDA）、卡爾曼濾波（Kalman Filtering）等經(jīng)典機器學(xué)習(xí)算法已應(yīng)用到腦電解碼當(dāng)中［8］。深度學(xué)習(xí)算法因其模型可以在多種場合下應(yīng)用，且不需要基于先驗經(jīng)驗知識的特征提取便能得到較好識別效果，從而獲得較多研究人員的青睞?；谛〔ǚ治龅姆椒商岣邊?shù)準確率，Neethu Robinson［9］通過腦電信號解碼獲得手部運動軌跡。精細回歸樹算法簡單實用，在腦電信號解碼過程中可以發(fā)揮較好的效果，是流行的腦電信號解碼方法［10］。

根據(jù)實驗獲得樣本特征，選擇合適且有效的數(shù)據(jù)解碼器，同時提出一種簡單可靠的腦電信號預(yù)處理和特征提取方法和腦電信號解碼策略，建立腦電信號和步態(tài)間的映射關(guān)系，是至關(guān)重要的。本文提出采用低通巴特沃斯濾波器對信號進行提取，并采用短時傅里葉變換和重采樣的連續(xù)小波變換方法對信號進行特征提取。提取后的腦電信號特征經(jīng)過精細回歸樹模型與步態(tài)信號建立聯(lián)系，實現(xiàn)了基于腦電信號的大鼠自由運動狀態(tài)下的步態(tài)解碼。

2 信號采集

為采集自由活動下大鼠的腦電和步態(tài)信號，搭建了專用實驗臺。實驗裝置由小動物走步機、支架、腦電采集裝置、攝像機、遮光筒等組成，如圖1所示。實驗過程中，大鼠通過“馬甲”固定于支架上，上肢懸空不負重，后肢負重在走步機履帶上且可以自由行走。大鼠行走速度取決于走步機轉(zhuǎn)速與大鼠后肢走路熟練程度。走步機可以通過調(diào)節(jié)履帶轉(zhuǎn)速調(diào)整大鼠行走速度。在大鼠行進過程中可通過此掛臂保持前臂懸空，只依靠后腿行進。遮光筒用于遮擋大鼠眼部周圍光線，減少外部環(huán)境干擾，使其行走時更加專注，同時降低腦電信號噪聲干擾。

圖1 實驗臺布置與信號采集示意圖

電極植入過程中，準確定位嚙齒動物模型腦區(qū)皮質(zhì)中后肢功能區(qū)域是重要的，這將決定信號采集的正確與否。本文采用Frost等人確定的后肢運動腦區(qū)對應(yīng)區(qū)域，定位其在前腦后部和顱骨縱向縫合中線的外側(cè)位置，從而確定大鼠大腦皮質(zhì)后肢運動表現(xiàn)的位置［11］。所有手術(shù)均在2～3%異氟醚富氧全麻下采用無菌技術(shù)完成。通過捏腳趾和控制呼吸頻率來監(jiān)測麻醉深度。動物被放在一個立體定位的框架上，手術(shù)過程中保持麻醉狀態(tài)并保證呼吸暢通。手術(shù)的第一步是在頭皮正中切開一道切口，從頭皮上撥開所有組織暴露出頭骨。然后，確定Bregma點和Lambda點，并標記所需的開顱位置。然后，在Lambda點的后部放置一枚螺絲釘以連接接地線，并在頭皮上額外安裝6枚螺絲釘以穩(wěn)定頭部的acrylic cement。微絲電極利用直徑在25 μm的商用金屬絲（鉑（70%）銥（30%）合金）制備。在M1的后肢區(qū)植入線間距為5 0 0 μm(1.5×1.5）的微絲陣列，電極為4 4排列，共16通道。在植入過程中將采集設(shè)備與電極相連，調(diào)整電極插入深度，并持續(xù)觀測數(shù)據(jù)記錄狀態(tài)，直到信號信噪比相對較高時停止。該陣列的中心被植入在Bregma后方2.6 mm，中線外側(cè)2 mm，硬腦膜表面下約1.5 mm的位置，覆蓋后肢區(qū)的所有區(qū)域。最后用牙科丙烯酸樹脂封閉暴露的開顱口和開口區(qū)，并做好術(shù)后護理。

實驗環(huán)境為一個相對黑暗狀態(tài)，大鼠腿部關(guān)鍵關(guān)節(jié)外側(cè)貼有反光球，走步機兩側(cè)各固定有用于追蹤大鼠步態(tài)的2臺3D攝像機。攝像機鏡頭邊緣附有補光燈，用于加強大鼠腿部關(guān)鍵點的反光效果，增強大鼠和腿部點的對比度，便于追光定位。步態(tài)信號通過Plexon CinePlex(Plexon Inc,Dallas,USA）步態(tài)采集系統(tǒng)進行采集，采集頻率為80 Hz。腦電信號采集系統(tǒng)（Plexon Inc,Dallas,USA）與步態(tài)信號采集裝置同步工作，信號通道數(shù)為16個，采樣頻率為40000 Hz。最終，將采集的步態(tài)數(shù)據(jù)與腦電數(shù)據(jù)存儲于PC中。

3 信號采集

3.1 腦電信號預(yù)處理

初步濾波后的大鼠腦電信號仍包含較多噪聲干擾，這些干擾信號在采集系統(tǒng)中16通道的腦電信號中同時含有，這種干擾項會降低腦電信號的信噪比。為減少此類噪聲的干擾，可通過均值降噪對信號進行消噪處理［12］。

16通道的腦電圖都包含一個參考信號，從原始數(shù)據(jù)中減去參考信號時，這個差值可以更清楚地反應(yīng)大腦活動的特征。使用信號幅值平均值作為參考信號Common Average Referencing(CAR）的方式可有效降低各個通道間信號噪聲。該方法一般用在多通道信號中，利用均值法處理腦電信號時，假設(shè)有N個通道的腦電信號，利用在某時刻N個通道的信號值總和計算平均值，通過減法將N個通道的數(shù)值減去平均值，獲得該通道的實際數(shù)值，即所需信號，其原理如公式（1）所示。

式中，N為信號通道數(shù)；yn為第n通道的腦電信號；y為原始腦電信號；yfiltered為均值濾波后的信號。本文中，16通道腦電信號均滿足信號質(zhì)量要求，對每個通道的信號做均值濾波處理，可獲取16通道的腦電信號，即n=16。

處理過程中，每個通道的信號均選用巴特沃斯濾波器進行濾波，以通道1信號為例，結(jié)果如圖2所示。橫坐標為時間，縱坐標為信號幅值。由圖中可看出，經(jīng)過均值濾波預(yù)處理之后，原始信號中高頻噪聲信號被消除，低頻信號得以突出，信號中的噪聲被進一步弱化。

圖2 均值濾波后第一通道信號

實驗過程中采集的腦電信號因采集頻率過高會在特征提取時產(chǎn)生過擬合。為采集有效的高質(zhì)量信號特征，需對濾波降噪后的腦電信號進行平滑處理。本文中，選用Savitzky-Golay濾波器（S-G濾波器）對均值濾波后的信號進行平滑處理。處理過程中，采用3001長度的3階窗函數(shù)S-G濾波器對信號進行濾波處理，結(jié)果如圖3所示。與未濾波信號相比，處理過的信號保留了足夠的原始信號的信息量，同時降低了信號波動和起伏程度，一定程度上提高了信噪比，有利于信號特征提取。

圖3 第一通道LFP信號通過S-G濾波平滑處理后結(jié)果

3.2 步態(tài)信號預(yù)處理

實驗過程中，大鼠腿部的標記點共有臀、膝、踝和足4個關(guān)節(jié)點位。步態(tài)信號通過信號采集系統(tǒng)獲取大鼠腿部關(guān)節(jié)點位的三維坐標。通過每個時間點獲取的坐標點位，每兩個相鄰關(guān)節(jié)點間可構(gòu)成一組向量。兩個向量之間的夾角可以反映某一關(guān)節(jié)隨時間變化的角度，本文以膝關(guān)節(jié)為參考點。向量的構(gòu)成方式如公式（2）所示：

式中，uk-h為膝部點位到臀部點位的向量；xk,yk,zk為膝關(guān)節(jié)三維坐標值；xh,yh,zh為臀關(guān)節(jié)三維坐標值。

基于關(guān)節(jié)點位向量可得出膝部點位到踝部點位、臀部點位的向量組?；谙嗤瑫r刻的向量組，可計算出以膝關(guān)節(jié)為頂點的夾角變化值，從而獲得大鼠膝關(guān)節(jié)沿Z方向的速度變化量。通過計算可獲得大鼠后肢Z軸準確位置，結(jié)果如圖4所示。

圖4 大鼠腿部膝關(guān)節(jié)速度隨時間變化

4 步態(tài)解碼方法

本文基于回歸樹線性模型進行建模，模型如公式（3）所示：

已知訓(xùn)練集D=(x1,y1),(x2,y2)…(xN,yN),其中X為輸入變量，f(x）為輸出變量。給定一組輸入空間劃分為R1,R2,R3…RM的數(shù)據(jù)集，每一個劃分的空間Rm對應(yīng)一個輸出值cm。在對空間劃分后回歸樹模型的優(yōu)劣程度可以用平方誤差來描述，如公式（4）所示：

每組數(shù)據(jù)輸出值的最優(yōu)解可描述為整組數(shù)據(jù)所有點的平均值，如公式（5）所示：

已有空間劃分情況下回歸樹的最優(yōu)生成可由公式（5）進行描述。一般情況下，在空間劃分過程中可采用以下方法：1）在所有輸入變量中找到第j個輸入變量xj，和它對應(yīng)的數(shù)值s作為初始的切分變量和切分點，并通過此輸入變量將輸入空間劃分成如式（6）和（7）所示的兩個數(shù)值空間：

基于空間劃分，可以獲得公式（8）的優(yōu)化函數(shù)用于尋找最優(yōu)切分點：

式中，c1,c2分別為初始子數(shù)據(jù)集R1和R2上的輸出值，min為最小值。

在建立回歸樹模型前將腦電和步態(tài)數(shù)據(jù)分成兩組數(shù)據(jù)集（各為50%），腦電特征作為預(yù)測數(shù)據(jù)集，步態(tài)數(shù)據(jù)運動學(xué)特征作為響應(yīng)數(shù)據(jù)集。在建立模型前，訓(xùn)練集數(shù)據(jù)特征采用主成分分析法（Principal Component Analysis,PCA）進行處理。主成分分析法是一種常用的數(shù)據(jù)分析方法［13,14］，這種方法的原理是對數(shù)據(jù)做線性變換，變換結(jié)果是一組各維度線性無關(guān)的數(shù)據(jù)表示形式，最終可以降低數(shù)據(jù)維度，提高訓(xùn)練成熟度。

5 實驗結(jié)果與討論

本文通過對局部場電位信號進行低通濾波處理，選取信噪比最高、特征最明顯的結(jié)果。之后將LFP信號按照不同頻率劃分為6個信號區(qū)間，采用不同的提取方式進行信號特征提取，基于提取的不同信號特征，利用回歸樹模型進行大鼠步態(tài)解碼預(yù)測，最終得出預(yù)測結(jié)果。

首先對16通道的信號質(zhì)量進行評估。在信號采集的16通道中，選擇單個通道數(shù)據(jù)進行檢驗。在信號比較中，從所列出的濾波方法、特征提取方法和不同的步態(tài)特征中選擇一種作為對所有通道數(shù)據(jù)的處理方法，分別對每一通道數(shù)據(jù)進行解碼處理。

為選取最優(yōu)信號，引入皮爾遜相關(guān)系數(shù)作為對于結(jié)果優(yōu)劣的判定條件。相關(guān)系數(shù)是通過計算兩組變量的協(xié)方差并歸一到方差平方根上，得到對于兩組變量相似性判斷的量化分析，定義如（9）式所示：

式中，Cov(X,Y）為變量X和Y的協(xié)方差；Var［X］和Var［Y］分別為變量X和Y的方差。

基于公式（9），對16通道單通道進行相關(guān)性計算，結(jié)果如表1所示，16個通道與步態(tài)信號均有在0.2-0.4之間的相關(guān)系數(shù)，體現(xiàn)出一定的正相關(guān)性。

表1 信號解碼膝關(guān)節(jié)位置坐標的皮爾遜相關(guān)系數(shù)

選取相關(guān)系數(shù)最高（效果最優(yōu)）的第四通道信號，基于不同的頻率區(qū)間將LFP信號進行劃分，對每一段信號分別進行均值處理和S-G降噪，然后通過小波變換對信號進行進一步的特征提取，通過回歸樹方法探究腦電信號和大鼠步態(tài)之間的映射關(guān)系。

如表2所示，為腦電信號不同信號區(qū)間段與大鼠步態(tài)信息相關(guān)度，從表中可得腦電信號低頻區(qū)在解碼效果上有著更好的效果。隨著信號頻率的增加，基于腦電信號的步態(tài)解碼效果逐步降低，γ頻段解碼效果最差。

表2 不同LFP信號區(qū)間段對大鼠步態(tài)解碼結(jié)果

選擇以300 Hz為截止頻率的低通巴特沃斯濾波器提取局部場電位信號。利用16通道信號所定義的均值濾波器進行降噪處理，并通過S-G濾波器平滑處理，分別經(jīng)過短時傅里葉變換、小波變換和希爾伯特-黃變換（HHT）三種方法對信號進行特征提取。1）濾波后的信號經(jīng)過窗長度為2000，窗交疊長度為1000的窗函數(shù)并結(jié)合窗內(nèi)傅里葉變換提取出短時傅里葉變換后的信號特征，將每個通道的501×7231個特征轉(zhuǎn)換為分貝值。將該分貝值作為信號特征進行擬合和解碼；2）對同一信號采用小波變換提取特征，從基于廣義“Morse”小波進行分解后的信號中提取出189×3616408維度的特征；3）利用希爾伯特-黃變換進行特征提取，單通道信號被分解成10個本征模函數(shù)，通過10個本征模函數(shù)對原始信號進行分解，生成101×3616408維度的信號特征。

解碼過程中，利用建立的精細回歸樹模型為基礎(chǔ)，通過50%數(shù)據(jù)的訓(xùn)練集搭建映射模型。為建立訓(xùn)練集與步態(tài)數(shù)據(jù)的一一對應(yīng)關(guān)系，需統(tǒng)一數(shù)據(jù)維度。短時傅里葉變換提供的特征維度與步態(tài)數(shù)據(jù)一致，可直接進行模型訓(xùn)練。小波變換提取的信號特征與步態(tài)數(shù)據(jù)不一致，分別采用1）插值法對步態(tài)數(shù)據(jù)進行補足的方式和2）腦電信號重采樣降低腦電信號維度的方式，使步態(tài)信號與腦電信號維度一致；由于小波變換中重采樣的方式可以提高測試集相關(guān)度，希爾波特-黃變換僅利用重采樣的處理方式降低腦電信號采樣頻率，使之與步態(tài)信號維度相同。表3所示為三種不同濾波方式在大鼠LFP信號對步態(tài)特征解碼中的效果。

表3 三種LFP信號特征提取方式解碼效果對比

通過對比各種信號特征提取方式，基于腦電信號重采樣的小波變換方法提取的腦電信號特征，通過精細回歸樹模型訓(xùn)練之后，在測試集中得到最好的結(jié)果，如圖5所示。預(yù)測所得結(jié)果和真實信號之間有著較高相關(guān)度，呈現(xiàn)正相關(guān)關(guān)系。但在信號幅值維度仍有一定差異，其主要原因為：小波變換作為一種精度較高的特征提取方式，在模型訓(xùn)練過程中容易產(chǎn)生過擬合現(xiàn)象，導(dǎo)致訓(xùn)練結(jié)果在趨勢上一致，在幅值上與步態(tài)信號產(chǎn)生差異。

圖5 LFP信號經(jīng)CWT特征提取解碼訓(xùn)練集結(jié)果

6 結(jié)論

為建立大鼠腦電信號和步態(tài)信號之間映射關(guān)系，解碼預(yù)測大鼠步態(tài)，提出了基于局部場電位信號的特征提取和步態(tài)解碼方法。在解碼方法探索過程中，發(fā)現(xiàn)采用300 Hz為截止頻率的低通巴特沃斯濾波器對局部場電位信號進行提取，并應(yīng)用短時傅里葉變換和重采樣的小波變換方法對信號進行特征提取的方式效果較好，測試集信號解碼相關(guān)度可達0.5081。提取后的腦電信號特征經(jīng)過精細回歸樹方法可以更好的與步態(tài)信號建立聯(lián)系。同時，探討了大鼠局部場電位信號與步態(tài)特征相關(guān)度的問題，發(fā)現(xiàn)頻率越低的信號（δ）對步態(tài)信號的解碼預(yù)測效果越優(yōu)，頻率越高的信號（high γ）對步態(tài)信號解碼預(yù)測效果越差。在接下來的研究中，更應(yīng)該關(guān)注低頻信號的采集和處理，這會給基于腦電信號為輸入的控制方法帶來幫助，促進腦-機接口技術(shù)的工程應(yīng)用。