基于延時等概率符號化傳遞熵分析的腦肌耦合雙向神經信息傳遞規律研究

張凱,徐光華,2,李文平,江開元,田沛源,鄭小偉,韓丞丞,張四聰

(1. 西安交通大學機械工程學院,710049,西安; 2. 西安交通大學機械制造系統工程國家重點實驗室,710054,西安)

隨著我國老齡化程度的持續加深,腦卒中疾病和意外腦損傷引發運動功能障礙的患者數量大幅增加[1],而手作為人體的精細感知和運動交互器官,在大腦皮層中面積占比大、功能復雜,因此手部運動功能極易受到損傷且難以修復[2]。神經接口作為一種融合生命科學與工程技術的前沿技術,能夠通過有效提取患者神經信號中包含的運動意圖。利用信號處理將其數字化為外部設備的控制指令,以輔助患者實現對外部設備的控制,從而實現主動運動意念參與的康復訓練,促進其手部運動功能恢復,因而得到了廣泛的研究和應用[3-5]。其中,腦電信號(EEG) 能夠從中樞層面反映患者殘存神經的運動意圖[6],而肌電信號(EMG) 能夠體現肌肉的活動狀態,體現了四周神經對于運動任務決策的神經響應[7]。同步耦合EEG與EMG不僅能夠更加全面地反映大腦對于肢體運動控制神經環路的工作機制,以有效指導神經接口界面的設計,從而幫助有一定自主性運動能力的患者進行手部精細動作的康復訓練[8]。同時,可以利用兩種信號的同步耦合關系來分析運動控制通路的功能狀態,為腦卒中患者的運動機能進行客觀量化的評價[9]。

然而,現有的研究只關注腦肌電信號模式識別的問題,忽略了對腦肌耦合的內在作用規律和雙向傳遞機制的研究,從而導致現有腦肌耦合神經接口存在著可解釋性差,控制準確率低等問題。腦肌電耦合通常指的是大腦皮層的運動控制信息會引起相關肌肉組織的同步振蕩活動[10-11],且該種振蕩存在雙向性,研究腦肌電耦合分析方法,不僅能夠有效揭示運動控制環路的神經信息內在作用機制,且能為高魯棒性的多層次腦肌耦合神經接口模型提供理論支撐[12]。現有應用于腦肌電信號的耦合分析方法從原理上可以分為線性方法和非線性方法[13-14]。其中,以相干性分析、格蘭杰因果分析為主的線性分析方法雖然具有計算簡便、可靠性高的優點,但不能滿足對腦肌電信號耦合的非線性特征的提取與辨識[15-17]。以互信息為主的非線性方法雖然可以對腦肌電的非線性過程進行刻畫,但是不能衡量系統狀態變化所帶來的動態信息變化。傳遞熵[18]通過計算轉移概率衡量系統間的動態信息傳遞,能夠分析時間序列間的非線性信息傳遞關系,是衡量兩個時間序列信息定向傳遞的指標,為腦肌電信號的耦合分析提供重要研究手段[19-20],然而,傳遞熵計算前需通過相空間提取和還原信號的動態特性,因此要求信號的長度足夠長,且傳遞熵的計算復雜度卻較高,極大限制了其應用。

基于此背景,本文面向腦肌耦合雙向神經信息傳遞規律開展研究,結合腦肌電信息傳遞的特征,將時間序列符號化的方法與傳遞熵計算相結合,提出了延時等概率符號化傳遞熵方法,以實現在粗粒化信號序列降低計算復雜度的同時保留腦肌電信號的動態特性,避免符號化過程中未考慮信號中出現幅值較大的突變噪聲會造成區間劃分不恰當,導致少數符號代替大部分原始觀測值的問題,以防損失大量原始時間序列中的信息。同時,基于所提方法開展在線實驗驗證,以揭示從腦肌電層面上的大腦皮層和運動肌肉間的雙向信息傳遞內在作用機制,為強魯棒性神經接口設計提供有力的理論支撐。

1 腦肌融合的神經傳導功能延時等概率符號化傳遞熵分析

傳遞熵通過衡量系統之間的信息交換來獲得系統的結構化信息,能夠獲得系統中的動態信息和信息交換的方向性,該方法能夠衡量兩個時間序列間的雙向非線性耦合特征,是研究腦電與肌電間信息交互的理想方法。為了避免傳遞熵相空間提取和還原過程的計算復雜度,研究人員將時間序列符號化的方法與傳遞熵計算相結合[21],一種典型的變尺度符號化的方法如下

R=

(1)

式中:R為符號化后的序列;x(i)為原始信號序列;k為符號化的數量;δ為遞增常量,其取值為(max(x)-min(x))/k。

基于該思路,研究人員提出變尺度符號化傳遞熵方法對腦肌電信號進行同步分析[22],其數學表達式為

(2)

式中:xn+τ為時間延遲τ時的信號X符號化觀測值;xn為信號X的符號化觀測值;yn為信號Y的符號化觀測值。

符號化傳遞熵能夠保留系統的動態特性并簡化傳遞熵的計算,由于其在符號化過程中未考慮信號中出現幅值較大的突變噪聲會造成區間劃分不恰當、導致少數符號代替大部分原始觀測值的問題[23],損失了大量原始時間序列中的信息。同時,腦肌電信號間的耦合關系并不僅是簡單直接的信息傳遞,由于大腦皮層到手部肌肉之間通過神經通路連接,而神經元信息傳遞存在延遲,可以推斷大腦與肌肉間的控制信息傳遞和感知反饋間必然也存在一定延遲[24-25]。傳統的傳遞熵分析方法僅考慮兩個時間序列間信息定向傳遞的方向,并未考慮兩個時間序列之間存在時間延遲的情況,所以直接使用傳統的傳遞熵方法對同步采集的腦肌電進行同步性分析并不能準確表達這種傳遞延遲。

因此,本文提出改進符號化傳遞熵計算方法,主要從符號化方法和傳遞熵的計算兩個方面進行改進。在符號化方法基礎上提出等概率符號化方法,在符號化區間劃分時盡可能保留信號的原始動態信息;在傳遞熵計算基礎上提出延時傳遞熵計算方法,在腦肌耦合分析時加入對腦肌電信號的信息傳輸延遲的考慮,算法流程如圖1所示。

圖1 基于延時等概率符號化傳遞熵的腦肌電信號耦合分析方法 Fig.1 Brain-muscle coupling analysis based on delay equal probability-symbolized transfer entropy

1.1 等概率符號化方法

本文所提等概率符號化方法,使符號化時間序列中各符號的出現概率相同,并可以根據需要動態調整符號數量。從而更好地保留原始時間序列中的信息,具體算法思路如下:對原始時間序列X(n)={x(1),x(2),…,x(n)}由小到大排序,排序后的序列為V(n)={v(1),v(2),…,v(n)};設定一個可變的尺度參數k,k表示需要劃分的符號數,則每個符號中對應的采樣點數j的計算式為

(3)

式中:round(·)為四舍五入運算。

取區間劃分參數

d(k)=v(kj)

(4)

式中:k為1,2,…,k-1。

然后,對原始時間序列X(n)進行符號化,符號化過程的函數形式為

(5)

由于n?k,因此不必考慮邊界點的歸屬和邊界區間極少數據點數差異,可以認為通過等概率劃分區間的方法,能夠保證每個區間內的數據點數近似相同。

從符號化時間序列包含的信息量角度進行考慮,符號時間序列中包含的信息量可以由香農熵表示,其計算式為

(6)

式中n為符號數量。可知當各個符號的出現概率p(xi)相同時,系統的不確定性和信息量最大。

傳統變尺度符號化方法對生理電信號進行符號化可能導致少數幾個符號代替大量的觀測值,因而其符號化時間序列中的不確定性較低,包含的信息量少,而等概率符號化方法使各個符號的出現概率相同,能夠包含原始時間序列中更多的動態信息,同時也具有傳統符號化方法降低信號中動力學噪聲和測量噪聲影響的優點。

1.2 延時等概率符號化傳遞熵

式(2)描述的方法僅對被預測樣本x(n+τ)添加延時參數,這種方式能夠一定程度上分析兩個時間序列的延遲,但它忽略了序列自身的信息傳遞,當時間序列間的信息傳遞延遲遠小于時間序列自身內的信息傳遞延遲時,這種改進傳遞熵方法并不能去除序列自身的信息傳遞影響,因而得不到準確的序列間信息傳遞結果。

本文所提延時等概率符號化傳遞熵(ETST),對于符號化時間序X(n)={x(1),x(2),…,x(n)}和Y(n)={y(1),y(2),…,y(n)},定義ETST(Y→X)為

ETST(Y→X)=

(7)

式中τ為序列Y(n)到X(n)的信息傳遞延遲。

按照傳遞熵的定義,本文定義的ETST(Y→X)的物理意義是已知信號Y和延遲τ后的信號X,去除信號X自身的信息量,信號Y能夠對延遲為τ+1的信號X的預測起到較大的作用。

2 腦肌電信號耦合分析實驗流程與數據采集

為了驗證所提方法的有效性,面向手部功能康復訓練中的常用手勢——手部抓握運動,開展手部運動過程中的EEG和sEMG信號的實時采集與分析實驗。實驗共招募了10名無精神類疾病的健康被試,7名男性和3名女性,均為右利手,平均年齡為24歲。所有被試者在實驗前均已閱讀實驗說明并簽署了知情同意書。每次實驗持續時間為12～17 s,具體流程如下:實驗開始,在t=0 s時,屏幕中央以文字提示被試者實驗即將開始,被試者得到提示后注視屏幕;以屏幕顯示靜止的左手和右手,屏幕中心顯示注視十字,被試者注視十字并保持靜息狀態,為排除被試自身節律或時間期待等因素的影響,靜息狀態持續時間隨機地設置為5～10 s;范式隨機提示左側或右側手部的運動想象任務,屏幕中心十字變為指向左或右的箭頭,根據范式提示,被試使用同側的手與范式進行持續5 s的同步抓握運動,且將抓握力保持在低于20%最大自主收縮的水平;單次任務結束,進入下一次的實驗準備階段。

本實驗使用gUSBamp設備同步采集腦電和肌電信號,腦電信號采集選取覆蓋運動皮層的腦電通道,具體包括F3、FZ、F4、FC3、FCZ、FC4、C3、CZ、C4、P3、PZ、P4;肌電信號采集兩側手臂的指淺屈肌(FDS)和肱肌(BM)處的肌電信號。

實驗要求被試者集中精力想象手部活動,實驗要求在安靜、微暗的環境下進行,避免外界噪聲對被試者造成干擾,實驗現場如圖2所示。最后獲得被試者每側手部以較小力度抓握時的腦肌電數據,選取每名被試者每側手部各20個試次,用來研究手部運動過程中大腦運動皮層和運動肌肉間的耦合關系。

圖2 腦肌電同步采集實驗環境示意圖Fig.2 Schematic diagram of experimental environment for synchronous EEG EMG acquisition

3 腦肌電傳遞延遲分析

對腦電信號進行1～100 Hz帶通濾波,以去除高頻噪聲和低頻運動偽跡的影響,然后使用延時等概率符號化傳遞熵對10名被試者的腦電信號和肌電信號進行分析,以獲得被試者的腦肌電間信息傳遞延遲時間。進行EEG至EMG的傳遞分析時,本文選使用等概率符號化方式對腦電信號和肌電信號分別進行符號化,設定參數k=10。

根據大腦對軀體的對側控制機制,選取右手指淺屈肌的肌電信號與C3通道腦電信號進行信息傳遞延遲分析,左手指淺屈肌的肌電信號與C4通道腦電信號進行信息傳遞延遲分析,現有研究表明皮層與肌肉的信息流延時大約在20～30 ms左右,因此本文選取1～50 ms延遲時長內的傳遞信息為分析對象,由此測量每個被試EEG至EMG信息流傳遞延遲時間,并對每個被試者的實驗數據多個試次進行疊加平均,從而完成進一步的觀察。其中,被試S1的兩側手部EEG→EMG傳遞熵隨延遲時間變化的分析結果,τ為最優延遲時間,如圖3、圖4所示。

圖3 左手EEG至EMG信息傳遞延遲分析Fig.3 The delay analysis of left hand EEG→EMG information transmission

圖4 右手EEG至EMG信息傳遞延遲分析Fig.4 The delay analysis of right hand EEG→EMG information transmission

可以看到,進行左側手部抓握時,在延遲時間為23.3 ms時出現傳遞熵峰值,在該時間延遲下腦電信號與肌電信號間存在最強的單向信息傳遞,可以認為是肌肉對大腦控制信號的響應延遲。在右手手部抓握時延遲時間為23.3 ms時出現傳遞熵峰值。

對EMG至EEG的信息傳遞延遲進行分析時,同樣采用C3、C4通道腦電信號與對側手部的指淺屈肌的肌電信號進行符號化傳遞熵的分析,結果如圖5、圖6所示。

圖6 右手EMG至EEG符號化傳遞熵分析Fig.6 The transfer entropy of right hand EMG to EEG information transmission

與EEG至EMG相似,被試S1的左手和右手抓握時EMG至EEG的信息傳遞延遲分別為29.2、30.8 ms。由于信息傳遞延遲可能存在一定波動,有的被試者可能出現峰值不顯著的情況,例如被試S1的右手EMG至EEG的傳遞熵峰值相對不顯著,但人體的神經傳導速度基本確定,在20～35 ms范圍內選擇峰值最可能符合實際信息傳遞延遲。記錄其他被試者的腦電與肌電的符號化傳遞熵分析結果,最優延遲時間的分析結果如表1所示。

表1 被試者S1～S10的最優延遲時間

由表1可知,在被試群體中,從EEG至EMG方向的最優延遲時間普遍短于從EMG至EEG,且EEG至EMG的延遲時長為21.7～32.5 ms,而EMG至EEG的時長約為24.2～33.2 ms。由于不同被試者之間存在個體差異,其最優延遲時間也存在差異,但該差異僅在毫秒之內,可以忽略不計,實驗得到最優延遲時間具有普適性。

4 多通道腦肌電符號化傳遞熵分析

為了分析手部運動過程中上臂肌肉與大腦運動功能區之間的耦合關系,本文對各被試者左、右手進行抓握時運動一側的指淺屈肌肌電信號與各導聯腦電信號和進行延時等概率符號化傳遞熵進行分析,延遲τ使用表1中確定的最優延遲時間,以被試S1為例,被試S1的延時等概率符號化傳遞熵分析結果如圖7所示。

(a)左手

(b)右手

由圖7可以發現,左側手部抓握時,C4通道的腦電信號與左手指淺屈肌的肌電信號存在最強的單向(EEG至EMG)信息傳遞,并且左手指淺屈肌的肌電信號與P4通道的腦電信號存在最強的單向(EMG至EEG)信息傳遞。右側手部抓握時,C3通道的腦電信號與右手指淺屈肌在兩個方向上的信息傳遞均是最強。對其他被試者的腦肌電符號化傳遞熵的分析也呈現類似的結果,例如對于被試S3,其左手和右手進行抓握時EEG至EMG方向信息傳遞最強的腦電通道分別為C4、FC3,這也恰好符合大腦的對側控制機制,實驗表明本文改進符號化傳遞熵能夠分析出腦電和肌電間的雙向信息傳遞關系,從而體現手部肌肉與大腦的運動皮層的雙向耦合關系。

同時可以觀察到腦電到肌電的延時等概率符號化傳遞熵顯著高于從肌電到腦電,從因果性的角度分析,可以認為腦電是因,肌電是果,主要的信息傳遞是由腦電到肌電。這與一般認知是相符的,當肢體運動發生時:一方面,大腦的運動皮層會發出信號,沿著運動傳導通路下行傳導至肢體的神經和肌肉;另一方面,本體感覺信息通過運動感知傳導通路上行傳導至大腦的運動感覺區域,大腦對感知到的運動信息加以分析后發出調節運動的指令,大腦在運動調控中起主導作用。

5 運動過程中的符號化傳遞熵變化分析

為了研究運動執行過程中的腦電和肌電的耦合強度變化情況,本文分析了靜息態到運動態的腦電及肌電信號符號化傳遞熵變化。具體步驟如下:選取多通道腦肌電耦合分析中得到的信息傳遞最強的腦電和肌電導聯的數據,進行去噪和1～100 Hz的帶通濾波處理;對單個試次的腦肌電數據進行分析,從單個試次靜息態開始,以1 s的時窗,以0.1 s間隔滑窗截取信號,分別對每段數據計算EEG至EMG和EMG至EEG的在表1中確定的最優延遲時間下的延時等概率符號化傳遞熵,構成符號化傳遞熵隨時間的變化曲線;分別對左右手的20個試次的符號化傳遞熵隨時間的變化曲線進行疊加,計算均值和方差,得到每個被試者左手和右手在抓握時的符號化傳遞熵變化曲線。對于被試S1和S2手部抓握過程的平均符號化傳遞熵分析結果如圖8所示,圖8中方差為實際計算結果的1/10。

(a)被試S1左手

(b)被試S2左手

(c)被試S1右手

(d)被試S2右手

圖8中黑色豎線表示屏幕提示被試者開始手部抓握運動的時刻,可得如下結論:在提示被試者開始手部抓握動作前,腦電與肌電信號的傳遞熵幅值均較低,EEG至EMG的符號化傳遞熵值僅略高于EMG至EEG,說明二者并未表現出直接的因果關系,且靜息狀態時大腦皮層與手部肌肉間的信息傳遞較少;在提示被試者開始手部動作1～2 s后,被試的EEG至EMG符號化傳遞熵值有顯著上升趨勢,在運動結束后上升趨勢逐漸消失;被試者EMG至EEG的符號化傳遞熵值有一定上升趨勢,但幅度較小,不如EEG至EMG的變化明顯,能夠說明在手部主動運動時,大腦皮層向肌肉發送調控指令,從而引起肌肉的活動,肌肉活動也會通過上行的通道向大腦皮層發送感知反饋信息,從而引起肌電到腦電方向的耦合強度上升;在兩側手部均使用較小抓握力度的情況下,各被試者右手EEG至EMG的符號化傳遞熵值要明顯高于左手,推測是由于本文被試者均為右利手,右手肌肉與大腦皮層的耦合性更強,一定程度上說明了本文符號化傳遞熵方法能夠較好地分析出腦肌耦合強度。

綜上所述,本文改進符號化傳遞熵方法可以分析大腦皮層和運動肌肉間的雙向信息傳遞,能夠體現腦肌耦合的時間延遲和耦合強度。經過分析,在人體手部以較小力度運動時時EEG至EMG方向的信息傳遞要高于EMG至EEG,并且不同方向上主動運動時的耦合強度均高于靜息狀態時。這可能能夠作為識別運動是由患者主動發起還是由于被試者的手部肌肉痙攣或其他癥狀引起的病態運動的依據,對于臨床康復訓練有一定的價值。

6 結 論

本文圍繞人體手部抓握運動過程中腦肌電神經信息的耦合規律開展研究,提出了延時等概率符號化傳遞熵分析方法,通過在線實驗同步采集手部運動任務中的EEG與EMG信號,進行了手部抓握過程中腦電信號和肌電信號間的雙向神經信息傳遞延遲檢測,開展了多通道的腦肌電耦合過程中的功能區定位與強度變化分析。通過延時等概率符號化傳遞熵對采集得到的EMG與EEG進行分析,可以得到如下結論。

(1)等概率符號化方法使各個符號的出現概率相同,能夠包含原始時間序列中更多的動態信息,可以更加精確的腦肌電信號的耦合特征,延時化的傳遞熵能夠有效能夠準確地檢測到兩個時間序列間信息傳遞,更能準確定位其信息傳遞延遲,在兩個序列間不存在信息傳遞時計算結果也保持了較好的魯棒性。

(2)通過對人體手部抓握運動的腦肌電信號進行分析,得到人體左右手的腦肌信息傳遞具有不對稱性,且該傳遞時延約為20～35 ms。左側手部抓握過程,C4通道的腦電信號與左手指淺屈肌的肌電信號存在最強的單向(EEG至EMG)信息傳遞,而右側手部抓握時,C3通道的腦電信號與右手指淺屈肌在兩個方向上的信息傳遞均是最強,且從腦電到肌電的傳遞熵峰值顯著高于從肌電到腦電的,指明了大腦在運動調控中起主導作用。同時,運動執行任務過程中,腦肌電耦合強度的變化規律表明了不同方向上主動運動耦合強度均高于運動靜息狀態。

綜上,本文不僅提出了全新的腦肌耦合分析技術方法并證明了該方法的有效性,同時揭示了抓握運動過程大腦皮層和手部肌肉的神經作用機制與雙向耦合關系,為新的康復方案和康復評價方法提供了理論依據,為神經接口技術的發展提供有力的算法支撐。