針對神經外科手術中用皮層電刺激(ECS)法確認病灶周邊運動功能區時存在操作時間長、可能損傷皮層以及誘發癲癇等缺陷,本研究提出一種基于皮質腦電(ECoG)分析的皮層運動功能區定位方法。該研究先對ECoG進行功率譜估計,分析Mu節律在運動前后的特征;再采用多分辨率小波分析方法提取特征并量化,觀察各導聯特征總量的分布,得到各導聯電極下皮層與運動功能區的相關程度。分析結果顯示,各電極下皮層與運動功能區的相關程度的分布與ECS的結果基本一致,且此定位方法的總操作時間大約為5 min。這表明基于ECoG的被動大腦皮層運動功能區定位方法,很可能成為術中ECS法的一種有效輔助定位手段。
引用本文: 符瓊琳, 姜濤, 黃岳山. 基于皮質腦電的大腦皮層運動功能區定位. 生物醫學工程學雜志, 2015, 32(4): 881-886. doi: 10.7507/1001-5515.20150157 復制
引言
由于人類大腦功能區的位置以及分布范圍存在個體差異性,因此在神經外科手術中,如何快速、準確、安全地確認病灶周邊的重要功能區,一直是外科醫生需要面臨的重要問題[1]。目前最常用的術中大腦皮層功能定位方法是皮層電刺激(electrocortical stimulation,ECS),且一直被視為大腦皮層功能定位方法中的“金標準”[2]。但ECS方法存在操作時間長、可能損傷大腦皮層以及誘發癲癇等問題[3-4],當患者為兒童時這些問題尤為突出;另外,ECS方法是一種抑制正常生理功能反應的方法,而非檢測方法[5]。為了改善直接電刺激法在定位運動功能區時存在的這些主要問題,并從檢測角度出發,本文提出一種基于皮質腦電(electrocorticography,ECoG)分析的大腦皮層運動功能區定位方法。
ECoG信號具有空間分辨率高、信噪比高,且長時間記錄時魯棒性好等特點[6],術中容易獲取。已有研究表明,ECoG信號在不同頻帶中包含多種有用的節律,并隱藏豐富的生理神經活動信息[7]。當對象進行不同活動時,大腦皮質的相應功能區會處于工作狀態,從而對應產生不同特點的ECoG,而且存在一定的頻域分布差異性;當人們肢體在真實運動或者想象肢體運動時,在來自運動皮層區域的ECoG信號中,Mu和beta節律會出現事件相關去同步化(event-related desynchronization,ERD)和事件相關同步化(event-related synchronization,ERS)[8]。鑒于小波分析的良好時頻分析特性[9],本研究采用小波分析方法從ECoG中提取出與運動相關的特征Mu節律,再量化其ERD現象,從而確定每一電極下大腦皮質與運動功能區的關系,完成運動功能區的定位。
1 數據采集
本研究所用ECoG數據來自廣州軍區廣州總醫院的兩名同意參加本項研究的癲癇病患者:對象A(男性)和對象B(女性),采集ECoG時所用電極陣列的電極個數分別為9和15。該數據的采集得到廣州軍區廣州總醫院倫理委員會的同意。
采集ECoG時,先在每位患者的大腦皮質上安放硬膜下電極陣列,患者重復做手掌張開-閉合動作4 s,然后休息4 s,重復進行2 min,腦電圖機同時采集并記錄ECoG信號,信號采樣頻率為250 Hz。
根據患者的運動、靜止狀況,分別從對象A的每一導ECoG信號中截取9個樣本,共有81個樣本;從對象B的每一導ECoG信號中截取9個樣本,共得135個樣本。
2 數據處理
2.1 特征選擇
當人在完成不同任務時,大腦皮質的相應功能區會處于工作狀態,從而產生不同特點的ECoG,而且存在一定的頻域分布差異性。根據此原理,對患者執行手掌張開/閉合任務時(即運動態)以及處在休息狀態(即靜息態)時的ECoG信號進行功率譜分析,結果如圖 1所示,其中圖 1(左)的樣本S1取自ECS法確認與運動相關的導聯,圖 1(右)的樣本S2取自ECS法確認與運動無關的導聯。從功率譜圖可知,圖 1(左)中,靜息狀態時的功率譜線在7~25 Hz之間比運動態的功率譜線高;這表明:與靜息狀態下相比,運動狀態下的Mu節律和beta節律在時域中的幅值會明顯減小,出現ERD現象,而且這種ERD現象在8~13 Hz頻帶(即Mu節律)表現尤為明顯。然而在圖 1(右)中,靜息狀態時的功率譜線與運動狀態時的功率譜線在7~25 Hz沒有明顯差異,所以,本文選擇Mu節律作為與運動相關的特征信號進行分析。

2.2 特征提取
ECoG是一種非平穩的隨機信號,近年涌現出多種分析非平穩隨機信號的方法,比如小波分析、獨立成分分量分析等。小波分析是一種時頻局部化分析的方法,在低頻部分具有較高的頻率分辨率,在高頻部分具有較高的時間分辨率[10],且它用可變的時頻窗口去分析信號的不同頻率成分:在大尺度下,可將信號的低頻信息全局特性表現出來,在小尺度下,可將高頻信息局部特性表現出來。所以小波變換對信號具有良好的自適應性。鑒于小波分析良好的時頻分析性能和自適應性,滿足既能分離出Mu節律,又能看得到Mu節律上ERD現象這一局部信息的需求,本研究采用多分辨率離散小波分析法提取與運動相關的特征信號Mu節律。
多分辨率離散小波分析法即Mallet算法,是Mallet提出多分辨率分析的概念時給出的一種小波分解的快速算法[11]。根據這一算法,若x(n)為連續腦電信號x(t)的離散采樣數據,x(n)=c0(k),則信號x(n)可以按式(1)分解:
$\begin{array}{l} {c_j}\left( k \right) = \sum\limits_n {h'} \left( {n - 2k} \right){c_{j - 1}}\left( n \right),\\ {d_j}\left( k \right) = \sum\limits_n {g'} \left( {n - 2k} \right){c_{j - 1}}\left( n \right), \end{array}$ |
其中cj和dj分別為j尺度上的逼近系數和細節系數;h′、g′為一對正交鏡像濾波器組;j為分解層數。
信號重構是信號分解的逆過程,其原理如式(2)所示,其中h、g分別是h′、g′的逆序。以第j尺度上的逼近系數cj和細節系數dj重構第j-1尺度上的逼近系數cj-1,直到c0,即得到原始信號。本研究選擇離散db3作為小波基,對ECoG信號進行4層小波分解,則各層子頻帶對應的頻率范圍如表 1所示,信號的采樣頻率fs為250 Hz。
${c_{j - 1}}\left( n \right) = \sum\limits_n {{h_{n - 2k}}{c_j}\left( k \right)} + \sum\limits_n {{g_{n - 2k}}{d_j}\left( k \right)} $ |
根據功率譜分析的結果,本研究選擇細節系數d4的單子頻帶重構信號Sd4作為分析信號,由表 1可知,其對應頻帶為7.812~15.625 Hz,包含與運動相關的特征信號Mu節律(即8~13 Hz頻帶);以該頻帶內的ERD現象作為運動區ECoG信號的特征,并按表達式(3)量化該特征[12]:
${\rm{ERD}} = \left( {{\rm{ER}} - {\rm{EA}}} \right){\rm{/ER,}}$ |
其中ER表示信號Sd4中運動時刻前800 ms內信號的能量;EA表示信號Sd4中運動時刻后800 ms內信號的能量。
統計每一導聯上樣本的特征量ERD,按表達式(4)計算各導聯上所有樣本的總特征量:
${K_n} = \sum\limits_i {{\rm{ER}}{{\rm{D}}_n}\left( i \right)} ,$ |
其中ERDn(i)表示來自第n導聯上第i個樣本的特征值,Kn表示第n導的總特征量,即所有樣本特征值的和,用于表征第n導聯電極下大腦皮層與運動功能區的相關程度:Kn越大,表示第n導聯電極下大腦皮層與運動功能區的相關程度越高,反之亦然。

3 結果與分析
對所有樣本進行4層離散db3多分辨率小波分解,重構細節系數d4對應的單子頻帶Sd4。圖 2是一段原始信號經4層離散db3多分辨小波分解后各單子頻帶的重構信號,其中Scj表示逼近系數cj的單子頻帶重構信號,Sdj表示細節系數dj的單子頻帶重構信號,如圖 2所示,細節系數d4的單子頻帶重構信號Sd4和其他單子頻帶重構信號相比,有明顯的ERD現象。

根據式(3)計算每一個樣本在Sd4上的特征值ERD,再按照式(4)計算每一導聯的總特征量Kn。各對象的總特征量分布情況如圖 3所示。

(a)對象A各導聯的總特征量分布;(b)對象B各導聯的總特征量分布
Figure3. Feature distribution(a) feature distribution of subject A; (b) feature distribution of subject B
圖 3(a)中,導聯1和導聯4對應電極下的大腦皮質與運動功能區的相關程度高,其次是導聯3和導聯2;圖 3(b)中,導聯9對應電極下的大腦皮質與運動功能區的相關程度最高,其次是導聯8和導聯11;而且各電極下大腦皮質與運動功能區的相關程度分布情況與各對象的直接電刺激結果基本一致。
這一結果在對各導聯的總特征量Kn之間的相關性分析中也有體現:如表 2所示,對象A的導聯1和導聯4相關程度極高,相關系數達到0.97,它們與導聯2和導聯3間的相關程度次之,而這兩導聯分別與其他導聯間的相關程度卻很低;同樣如表 3所示,對象B的導聯9與導聯8的相關系數為0.91,而與導聯11的相關系數為0.89,且這三個導聯與其他導聯間的相關系數均不高。這表明:運動功能區相關皮質產生的電信號的特征量Kn之間相關性強,而運動功能區相關的皮質和非運動功能區相關的皮質所產生的電信號的特征量Kn之間僅弱相關。


以上分析結果表明基于ECoG分析的大腦皮層功能區定位的分析方法具有可行性。
4 討論和結論
Mu節律在執行運動任務前后出現的ERD和ERS現象被認為代表著與運動丘腦皮層調制相關的突觸后電位的變化[13],同時也被視為感覺運動區的空閑節律。學者Adrian和Matthews在提出“空閑系統(idling system)”的概念時指出[14],大腦在既不接收也不處理感覺信息時處于空閑狀態,皮層的空閑程度可以反映出皮層是否參與了感覺運動信息的處理。當運動任務狀態和休息狀態交替出現時,感覺運動皮層上Mu節律這一空閑節律出現的ERD和ERS現象正反映皮層參與感覺運動信息的處理的情況。而目前在國內,把Mu節律及其ERD和ERS現象應用到神經外科臨床中的功能定位的研究還比較少。
基于ECoG分析的大腦皮層功能區定位方法,先采用小波分析法提取出與運動相關的特征信號Mu節律,然后統計每一導聯上所有樣本的特征量ERD的總和,從而反映每一電極下大腦皮質與運動功能區的關系。該方法從數據采集至基本確定電極下皮層與運動功能區的關系大約只需5 min,而研究的分析結果可以指導外科醫生用ECS法做進一步精確定位,減小了ECS法的盲目性,縮短了手術的時間,也減小了患者的痛苦。此外,對分析源(即ECoG信號)的檢測屬于被動檢測,避免了ECS法可能導致的大腦皮層損傷或誘發癲癇等情況,所以本研究所提出的方法能提高手術中定位運動功能區操作的安全性。
在所分析的兩份ECoG數據中,各電極下大腦皮質和運動功能區相關程度的分布與“金標準”的結果基本一致,初步研究結果表明基于ECoG分析的大腦皮層運動區功能定位的方法有可能成為術中ECS法的一種有效輔助定位手段。隨著更多關于ECoG應用于大腦功能定位研究的進行,我們期待基于ECoG分析的功能定位這一新興研究方法具有更大的臨床實際應用價值。
引言
由于人類大腦功能區的位置以及分布范圍存在個體差異性,因此在神經外科手術中,如何快速、準確、安全地確認病灶周邊的重要功能區,一直是外科醫生需要面臨的重要問題[1]。目前最常用的術中大腦皮層功能定位方法是皮層電刺激(electrocortical stimulation,ECS),且一直被視為大腦皮層功能定位方法中的“金標準”[2]。但ECS方法存在操作時間長、可能損傷大腦皮層以及誘發癲癇等問題[3-4],當患者為兒童時這些問題尤為突出;另外,ECS方法是一種抑制正常生理功能反應的方法,而非檢測方法[5]。為了改善直接電刺激法在定位運動功能區時存在的這些主要問題,并從檢測角度出發,本文提出一種基于皮質腦電(electrocorticography,ECoG)分析的大腦皮層運動功能區定位方法。
ECoG信號具有空間分辨率高、信噪比高,且長時間記錄時魯棒性好等特點[6],術中容易獲取。已有研究表明,ECoG信號在不同頻帶中包含多種有用的節律,并隱藏豐富的生理神經活動信息[7]。當對象進行不同活動時,大腦皮質的相應功能區會處于工作狀態,從而對應產生不同特點的ECoG,而且存在一定的頻域分布差異性;當人們肢體在真實運動或者想象肢體運動時,在來自運動皮層區域的ECoG信號中,Mu和beta節律會出現事件相關去同步化(event-related desynchronization,ERD)和事件相關同步化(event-related synchronization,ERS)[8]。鑒于小波分析的良好時頻分析特性[9],本研究采用小波分析方法從ECoG中提取出與運動相關的特征Mu節律,再量化其ERD現象,從而確定每一電極下大腦皮質與運動功能區的關系,完成運動功能區的定位。
1 數據采集
本研究所用ECoG數據來自廣州軍區廣州總醫院的兩名同意參加本項研究的癲癇病患者:對象A(男性)和對象B(女性),采集ECoG時所用電極陣列的電極個數分別為9和15。該數據的采集得到廣州軍區廣州總醫院倫理委員會的同意。
采集ECoG時,先在每位患者的大腦皮質上安放硬膜下電極陣列,患者重復做手掌張開-閉合動作4 s,然后休息4 s,重復進行2 min,腦電圖機同時采集并記錄ECoG信號,信號采樣頻率為250 Hz。
根據患者的運動、靜止狀況,分別從對象A的每一導ECoG信號中截取9個樣本,共有81個樣本;從對象B的每一導ECoG信號中截取9個樣本,共得135個樣本。
2 數據處理
2.1 特征選擇
當人在完成不同任務時,大腦皮質的相應功能區會處于工作狀態,從而產生不同特點的ECoG,而且存在一定的頻域分布差異性。根據此原理,對患者執行手掌張開/閉合任務時(即運動態)以及處在休息狀態(即靜息態)時的ECoG信號進行功率譜分析,結果如圖 1所示,其中圖 1(左)的樣本S1取自ECS法確認與運動相關的導聯,圖 1(右)的樣本S2取自ECS法確認與運動無關的導聯。從功率譜圖可知,圖 1(左)中,靜息狀態時的功率譜線在7~25 Hz之間比運動態的功率譜線高;這表明:與靜息狀態下相比,運動狀態下的Mu節律和beta節律在時域中的幅值會明顯減小,出現ERD現象,而且這種ERD現象在8~13 Hz頻帶(即Mu節律)表現尤為明顯。然而在圖 1(右)中,靜息狀態時的功率譜線與運動狀態時的功率譜線在7~25 Hz沒有明顯差異,所以,本文選擇Mu節律作為與運動相關的特征信號進行分析。

2.2 特征提取
ECoG是一種非平穩的隨機信號,近年涌現出多種分析非平穩隨機信號的方法,比如小波分析、獨立成分分量分析等。小波分析是一種時頻局部化分析的方法,在低頻部分具有較高的頻率分辨率,在高頻部分具有較高的時間分辨率[10],且它用可變的時頻窗口去分析信號的不同頻率成分:在大尺度下,可將信號的低頻信息全局特性表現出來,在小尺度下,可將高頻信息局部特性表現出來。所以小波變換對信號具有良好的自適應性。鑒于小波分析良好的時頻分析性能和自適應性,滿足既能分離出Mu節律,又能看得到Mu節律上ERD現象這一局部信息的需求,本研究采用多分辨率離散小波分析法提取與運動相關的特征信號Mu節律。
多分辨率離散小波分析法即Mallet算法,是Mallet提出多分辨率分析的概念時給出的一種小波分解的快速算法[11]。根據這一算法,若x(n)為連續腦電信號x(t)的離散采樣數據,x(n)=c0(k),則信號x(n)可以按式(1)分解:
$\begin{array}{l} {c_j}\left( k \right) = \sum\limits_n {h'} \left( {n - 2k} \right){c_{j - 1}}\left( n \right),\\ {d_j}\left( k \right) = \sum\limits_n {g'} \left( {n - 2k} \right){c_{j - 1}}\left( n \right), \end{array}$ |
其中cj和dj分別為j尺度上的逼近系數和細節系數;h′、g′為一對正交鏡像濾波器組;j為分解層數。
信號重構是信號分解的逆過程,其原理如式(2)所示,其中h、g分別是h′、g′的逆序。以第j尺度上的逼近系數cj和細節系數dj重構第j-1尺度上的逼近系數cj-1,直到c0,即得到原始信號。本研究選擇離散db3作為小波基,對ECoG信號進行4層小波分解,則各層子頻帶對應的頻率范圍如表 1所示,信號的采樣頻率fs為250 Hz。
${c_{j - 1}}\left( n \right) = \sum\limits_n {{h_{n - 2k}}{c_j}\left( k \right)} + \sum\limits_n {{g_{n - 2k}}{d_j}\left( k \right)} $ |
根據功率譜分析的結果,本研究選擇細節系數d4的單子頻帶重構信號Sd4作為分析信號,由表 1可知,其對應頻帶為7.812~15.625 Hz,包含與運動相關的特征信號Mu節律(即8~13 Hz頻帶);以該頻帶內的ERD現象作為運動區ECoG信號的特征,并按表達式(3)量化該特征[12]:
${\rm{ERD}} = \left( {{\rm{ER}} - {\rm{EA}}} \right){\rm{/ER,}}$ |
其中ER表示信號Sd4中運動時刻前800 ms內信號的能量;EA表示信號Sd4中運動時刻后800 ms內信號的能量。
統計每一導聯上樣本的特征量ERD,按表達式(4)計算各導聯上所有樣本的總特征量:
${K_n} = \sum\limits_i {{\rm{ER}}{{\rm{D}}_n}\left( i \right)} ,$ |
其中ERDn(i)表示來自第n導聯上第i個樣本的特征值,Kn表示第n導的總特征量,即所有樣本特征值的和,用于表征第n導聯電極下大腦皮層與運動功能區的相關程度:Kn越大,表示第n導聯電極下大腦皮層與運動功能區的相關程度越高,反之亦然。

3 結果與分析
對所有樣本進行4層離散db3多分辨率小波分解,重構細節系數d4對應的單子頻帶Sd4。圖 2是一段原始信號經4層離散db3多分辨小波分解后各單子頻帶的重構信號,其中Scj表示逼近系數cj的單子頻帶重構信號,Sdj表示細節系數dj的單子頻帶重構信號,如圖 2所示,細節系數d4的單子頻帶重構信號Sd4和其他單子頻帶重構信號相比,有明顯的ERD現象。

根據式(3)計算每一個樣本在Sd4上的特征值ERD,再按照式(4)計算每一導聯的總特征量Kn。各對象的總特征量分布情況如圖 3所示。

(a)對象A各導聯的總特征量分布;(b)對象B各導聯的總特征量分布
Figure3. Feature distribution(a) feature distribution of subject A; (b) feature distribution of subject B
圖 3(a)中,導聯1和導聯4對應電極下的大腦皮質與運動功能區的相關程度高,其次是導聯3和導聯2;圖 3(b)中,導聯9對應電極下的大腦皮質與運動功能區的相關程度最高,其次是導聯8和導聯11;而且各電極下大腦皮質與運動功能區的相關程度分布情況與各對象的直接電刺激結果基本一致。
這一結果在對各導聯的總特征量Kn之間的相關性分析中也有體現:如表 2所示,對象A的導聯1和導聯4相關程度極高,相關系數達到0.97,它們與導聯2和導聯3間的相關程度次之,而這兩導聯分別與其他導聯間的相關程度卻很低;同樣如表 3所示,對象B的導聯9與導聯8的相關系數為0.91,而與導聯11的相關系數為0.89,且這三個導聯與其他導聯間的相關系數均不高。這表明:運動功能區相關皮質產生的電信號的特征量Kn之間相關性強,而運動功能區相關的皮質和非運動功能區相關的皮質所產生的電信號的特征量Kn之間僅弱相關。


以上分析結果表明基于ECoG分析的大腦皮層功能區定位的分析方法具有可行性。
4 討論和結論
Mu節律在執行運動任務前后出現的ERD和ERS現象被認為代表著與運動丘腦皮層調制相關的突觸后電位的變化[13],同時也被視為感覺運動區的空閑節律。學者Adrian和Matthews在提出“空閑系統(idling system)”的概念時指出[14],大腦在既不接收也不處理感覺信息時處于空閑狀態,皮層的空閑程度可以反映出皮層是否參與了感覺運動信息的處理。當運動任務狀態和休息狀態交替出現時,感覺運動皮層上Mu節律這一空閑節律出現的ERD和ERS現象正反映皮層參與感覺運動信息的處理的情況。而目前在國內,把Mu節律及其ERD和ERS現象應用到神經外科臨床中的功能定位的研究還比較少。
基于ECoG分析的大腦皮層功能區定位方法,先采用小波分析法提取出與運動相關的特征信號Mu節律,然后統計每一導聯上所有樣本的特征量ERD的總和,從而反映每一電極下大腦皮質與運動功能區的關系。該方法從數據采集至基本確定電極下皮層與運動功能區的關系大約只需5 min,而研究的分析結果可以指導外科醫生用ECS法做進一步精確定位,減小了ECS法的盲目性,縮短了手術的時間,也減小了患者的痛苦。此外,對分析源(即ECoG信號)的檢測屬于被動檢測,避免了ECS法可能導致的大腦皮層損傷或誘發癲癇等情況,所以本研究所提出的方法能提高手術中定位運動功能區操作的安全性。
在所分析的兩份ECoG數據中,各電極下大腦皮質和運動功能區相關程度的分布與“金標準”的結果基本一致,初步研究結果表明基于ECoG分析的大腦皮層運動區功能定位的方法有可能成為術中ECS法的一種有效輔助定位手段。隨著更多關于ECoG應用于大腦功能定位研究的進行,我們期待基于ECoG分析的功能定位這一新興研究方法具有更大的臨床實際應用價值。