基于腦電信號的疼痛強度識別方法研究

李冬，金韜，馮智英，左路路，朱翔，劉偉明

1.浙江大學(xué)信息與電子工程學(xué)院，浙江杭州310027；2.浙江大學(xué)第一附屬醫(yī)院疼痛科，浙江杭州310027

前言

疼痛是集病理、情感為一體的不愉快的主觀感受［1］。由于人體對疼痛的感知涉及到大腦的多個區(qū)域，通過對大腦活動的研究解碼疼痛變得非常復(fù)雜。腦電圖（Electroencephalogram,EEG）反映人體大腦皮層電活動的平均值，包含大量的生理病理信息［2］。通過對腦電信號進(jìn)行特征提取，有望實現(xiàn)對疼痛的客觀評估，為疼痛的臨床治療提供更多信息。

國內(nèi)外學(xué)者針對疼痛腦電信號特征提取方法進(jìn)行大量研究，采用包括時域［3］、頻域［4-6］、時頻域［7-10］和非線性動力學(xué)［1,11］等研究方法。時域分析方法主要是分析腦電信號波形的幾何性質(zhì)，如腦電信號的幅值、方差以及均值等。頻域分析方法，主要是通過FFT變換獲得功率譜，將腦電信號各頻率段的功率、功率譜峰值頻率［4］等作為研究疼痛的潛在特征。由于腦電信號是一種非平穩(wěn)并且產(chǎn)生機(jī)理非常復(fù)雜的隨機(jī)信號，單純從時域或頻域提取特征僅能夠分析出大腦在不同疼痛狀態(tài)下腦激活程度的差異，不足以從多維度提取有效特征對疼痛等級進(jìn)行分類預(yù)測。小波變換屬于時頻分析方法，其能夠提供一個隨頻率改變的“時間-頻率”窗口，具有良好的時頻分辨率，能夠捕捉到疼痛隨時間變化的局部信息，已被成功應(yīng)用在疼痛研究領(lǐng)域［7-9］。Hadjileontiadis［7］通過讓健康受試者的右手深入冰水中獲得冷痛刺激，特征提取采用基于小波變換的雙譜分析，在不痛和疼痛兩分類情況下獲得84.12%的分類準(zhǔn)確率。Vijayakumar等［9］對受試者進(jìn)行熱痛刺激，特征提取采用獨立成分分析和小波變換相結(jié)合的方法，以隨機(jī)森林（Random Forest,RF）為分類器，在1～10級的疼痛識別中準(zhǔn)確率高達(dá)89.45%。這些方法對外部刺激誘發(fā)的疼痛識別取得較好的分類效果，但沒有考慮到對原始特征集進(jìn)行特征選擇。原始特征集通常維度較大，計算成本高，且具有高冗余性，進(jìn)一步優(yōu)化特征集可以提高分類性能。RF利用袋外數(shù)據(jù)能夠評估變量重要性，可以有效地用于特征選擇，已被成功應(yīng)用于生物信息學(xué)［12］、黑色素瘤診斷等領(lǐng)域［13］。從目前來看，在疼痛識別領(lǐng)域尚未發(fā)現(xiàn)采用RF篩選最優(yōu)特征子集。

本研究通過臨床采集疼痛患者腦電信號，彌補外部刺激不能充分模擬自然和臨床痛苦的缺陷。采用基于離散小波變換（Discrete Wavelet Transform,DWT）的疼痛腦電特征組合方法，并利用RF進(jìn)行特征選擇。對特征選擇前后的分類準(zhǔn)確率進(jìn)行對比，實驗結(jié)果表明進(jìn)行特征選擇后能以較高的準(zhǔn)確率識別出不同的疼痛水平。

1 基本原理

1.1 DWT

DWT被廣泛應(yīng)用于生物信號的時頻分析，腦電信號具有非平穩(wěn)特性，利用小波變換能夠捕捉大腦感知疼痛的細(xì)節(jié)信息［14］。在眾多小波家族中，Daubechies類小波在腦電信號的分解中具有很多優(yōu)點［15］。因此，在本研究中采用Daubechies 4（db4）小波對腦電信號進(jìn)行4層分解。EEG信號的具體分解過程如圖1所示，首先原始信號x[n]通過半帶高通濾波器g[n]和低通濾波器h[n]被分解為低頻和高頻部分，分別得到近似系數(shù)A1和細(xì)節(jié)系數(shù)D1，然后對低頻部分重復(fù)上述步驟直到達(dá)到所需分解的層數(shù)。其中，分解系數(shù)的表達(dá)式如式（1）和式（2）所示：

其中，ai[k]和di[k]分別表示信號在第i層分解中第k點變換得到的近似和細(xì)節(jié)系數(shù)。

1.2 特征提取

圖1 腦電信號4層小波分解Fig 1 Four-level wavelet decomposition of electroencephalogram(EEG)signals

與疼痛相關(guān)特征的提取對于疼痛的客觀評估至關(guān)重要。由于腦電信號在不同的疼痛狀態(tài)下，各個腦電節(jié)律的能量分布、信號復(fù)雜度會發(fā)生變化［11,16］。因此，我們根據(jù)EEG信號經(jīng)DWT變換后的分解系數(shù)，計算各子帶的能量占比、樣本熵以區(qū)分不同的腦電狀態(tài)。其中，能量占比計算過程參考文獻(xiàn)［17］中的方法。樣本熵值能夠衡量時間序列的復(fù)雜度，時間序列越復(fù)雜，其對應(yīng)的樣本熵值越大［18］。樣本熵算法的計算過程參考文獻(xiàn)［19］的方法。其中引入兩個重要的參數(shù)m和r，m代表嵌入維數(shù)，r代表時間序列標(biāo)準(zhǔn)差。研究表明，當(dāng)參數(shù)m取值為1或2、r取值為0.10～0.25時，樣本熵值有較好的統(tǒng)計效果［19］。本研究中，m取2，r取0.15。

相位同步能夠揭示大腦不同區(qū)域神經(jīng)電活動的關(guān)系，已被廣泛應(yīng)用于腦電信號分析［20］。為了提取在不同疼痛狀態(tài)下各通道腦電信號間的相位同步關(guān)系，我們提出一種基于小波分解的相位同步計算方法。對于8通道腦電信號，分別計算在不同分解水平和不同通道對之間的鎖相值，共得到5個特征向量，每個特征向量包含28個元素。小波鎖相值的計算過程如下：

（1）設(shè)信號x(t)和y(t)分別為通道x和通道y的腦電信號經(jīng)DWT第l(1～4)層的分解信號，對信號x(t)和y(t)進(jìn)行希爾伯特變換得到信號xh(t)和yh(t)：

（2）計算信號xh(t)和yh(t)在時間t處的瞬時相位，得到瞬時相位差：

（3）信號xh(t)和yh(t)的鎖相值（PLV）被定義為兩信號在所有時間點的相位差的平均值：

1.3 基于RF的特征選擇與分類

RF是被廣泛應(yīng)用于數(shù)據(jù)挖掘、生物信息處理等領(lǐng)域中解決高維度和非線性樣本的一種分類器［21］。本研究中，為了驗證特征提取方法的有效性和提高分類模型的預(yù)測精度，我們采用基于RF的Wrapper特征選擇算法。利用基于袋外數(shù)據(jù)分類準(zhǔn)確率的變量重要性度量方式［22］，對特征進(jìn)行重要性排序，然后采用序列后向選擇算法（Sequential Backward Selection,SBS），每次從特征集合中去掉一個重要性得分最小的特征，逐次進(jìn)行迭代，并計算分類準(zhǔn)確率，最終將分類準(zhǔn)確率最高的特征集合作為特征選擇結(jié)果。為了評估分類器的魯棒性，在每次迭代中，采用留一交叉驗證方法，每次只留下一個樣本作測試集，其他樣本被當(dāng)做訓(xùn)練集［23］。同時，為了測試分類器的性能，我們分別計算了分類準(zhǔn)確率、靈敏性和特異性。

2 數(shù)據(jù)采集與預(yù)處理

本研究中，腦電數(shù)據(jù)共包含3個數(shù)據(jù)集（無痛、輕痛和重痛），無痛組數(shù)據(jù)來自于7例健康受試者，輕痛和重痛數(shù)據(jù)集是從14例患有不同程度帶狀皰癥后遺神經(jīng)痛的患者中臨床采集到的腦電信號，根據(jù)患者對疼痛的自我描述，將數(shù)據(jù)歸類到輕痛組或重痛組。實驗設(shè)備采用8通道的腦電信號采集儀，采樣頻率為250 Hz，并自帶50 Hz陷波器和0～48 Hz的帶通濾波器，其電極安放位置根據(jù)國際10-20系統(tǒng)安放在C3、C4、P7、P8、O1、O2、FP1和FP2位置，參考電極固定在雙耳垂。在數(shù)據(jù)采集過程中，要求受試者閉眼放松，采集時間約3 min。對于明顯屬于干擾噪聲的信號，進(jìn)行手動剔除。最終，每組腦電數(shù)據(jù)共包含5 000個采樣點，用于后續(xù)特征提取。在實驗開始前，所有參與者均簽署知情同意書。

3 結(jié)果與分析

3.1 特征分析

首先利用db4小波對每位受試者的8通道腦電信號進(jìn)行4層小波分解，得到近似系數(shù)A4和細(xì)節(jié)系數(shù)D1～D4，具體分解結(jié)果如圖2所示。每層分解信號的子帶頻率范圍及其對應(yīng)的近似腦電節(jié)律如表1所示。基于每層分解的小波系數(shù)，我們計算小波統(tǒng)計特征值（最大值、最小值、平均值和標(biāo)準(zhǔn)差）、各子帶能量占比、樣本熵以及兩兩通道對之間的鎖相值。最終，對每組腦電信號，我們提取包含380個元素的特征向量。這些特征向量的具體描述如表2所示。

圖2 某患者原始EEG信號及DWT分解結(jié)果Fig 2 Original EEG signals from one patient and results of discrete wavelet transform(DWT)decomposition

表1 每層小波分解的子帶頻率范圍及其對應(yīng)腦電節(jié)律Tab 1 Sub-band frequency range of each wavelet decomposition level and the corresponding brain wave rhythms

表2 提取特征描述Tab 2 Descriptions of extracted features

3.2 分類結(jié)果

對提取的特征向量，我們采用基于RF的序列后向選擇（RF-SBS）算法進(jìn)行特征選擇和分類預(yù)測。疼痛強度的分類共包含3個兩分類問題（無痛和輕痛、輕痛和重痛、無痛和重痛）和1個三分類問題（無痛和輕痛和重痛）。同時，為了評估分類器性能，我們分別計算分類器的準(zhǔn)確度（Accuracy）、靈敏性（Sensitivity）和特異性（Specificity），各參數(shù)計算公式如下：

表3顯示由RF-SBS算法篩選后的特征作為RF分類器的輸入進(jìn)行訓(xùn)練和分類的預(yù)測結(jié)果。由表3可知，無痛和輕痛的分類準(zhǔn)確率為92.85%，無痛和重痛的分類準(zhǔn)確性為100%，輕痛和重痛的分類準(zhǔn)確率為92.85%，不同疼痛水平三分類的分類準(zhǔn)確率為80.95%。這表明所提出的方法可以作為疼痛水平識別的有效工具。

表3 RF-SBS特征選擇后的分類結(jié)果（%）Tab 3 Classification results using features selected by RF-SBS(%)

此外，我們將由RF-SBS算法選擇出的特征組合與所有特征分別作為分類器輸入時的分類結(jié)果做對比，結(jié)果如表4所示?？梢钥闯?，經(jīng)特征選擇后的分類效果明顯好于將所有特征作為分類器輸入的情況。這表明提取的腦電特征在表征和疼痛有關(guān)的信息之外，仍然具有較大的冗余性，不利于分類器模型的訓(xùn)練和預(yù)測。對高維特征進(jìn)行特征選擇后，疼痛預(yù)測的準(zhǔn)確度顯著提高。

表4 不同特征下的分類結(jié)果比較（%）Tab 4 Comparison of classification results based on different features(%)

4 結(jié)論

本研究提出一種新穎的方法用于識別帶狀皰疹后遺神經(jīng)痛患者的疼痛強度。基于DWT，對臨床腦電數(shù)據(jù)進(jìn)行特征提取，并使用RF-SBS算法對高維特征進(jìn)行特征選擇。這些被選特征被用于訓(xùn)練和預(yù)測疼痛水平，以21例臨床采樣數(shù)據(jù)作為訓(xùn)練和測試樣本，采用留一交叉驗證方式，達(dá)到較高的分類效果，平均準(zhǔn)確率為91.7%，其中對于無痛和重痛的分類情況達(dá)到100%的準(zhǔn)確率。這表明提取的特征能夠捕捉到由不同疼痛狀態(tài)引起的腦電模式的變化。另外，與將所有特征作為分類器輸入相比，經(jīng)RF-SBS篩選后的特征訓(xùn)練出的模型準(zhǔn)確度更高，說明高維腦電特征具有冗余性，篩選后的特征更能反映患者疼痛狀態(tài)。這種方法能夠作為一種有效的工具用來診斷疼痛水平，并為臨床疼痛的治療提供及時的信息。未來，可以考慮通過融合其他相關(guān)生物特征以進(jìn)一步提高分類準(zhǔn)確率，并能解釋潛在的疼痛反應(yīng)機(jī)制。