腦機接口系統的核心問題之一是信號分類。本文針對腦電信號的異構融合特征的分類問題提出了一種新方法：封裝式彈性網特征選擇和分類。首先，對預處理后的腦電（EEG）信號聯合應用時域統計、功率譜、共空間模式和自回歸模型方法提取高維異構融合特征。其次，采用封裝方式進行特征選擇：對訓練數據采用彈性網罰邏輯回歸擬合模型，通過坐標下降法估計模型參數，運用10倍交叉驗證選擇出最優特征子集。最后采用已訓練的最優模型對測試樣本進行分類。實驗中采用國際BCI競賽Ⅳ的EEG數據，結果表明，該方法適用于高維融合特征的最優特征子集選擇問題，對于EEG信號的識別不僅效果好、速度快，而且能夠選出與分類更相關的子集，獲得相對簡單的模型，平均測試正確率達到了81.78%。
引用本文： 李靜, 王金甲, 李慧. 融合腦電特征的彈性網特征選擇和分類. 生物醫學工程學雜志, 2016, 33(3): 413-419. doi: 10.7507/1001-5515.20160070 復制
引言
腦機接口(brain-computer interface,BCI)是大腦和計算機或電子設備之間傳遞信息的通訊系統，這個系統不依賴外圍神經和肌肉組成的大腦正常通路[1]。這種全新的人機交互方式能夠將大腦產生的信息轉換成對外部設備的控制命令，從而代替肢體和語言實現與外界的交流。這項技術不僅能夠為思維正常但有運動障礙的人提供幫助，還能作為特殊環境下正常人對外部設備的控制方式，甚至能夠提供一種全新娛樂體驗。
BCI系統的核心問題之一是腦電圖(electroencephalogram,EEG)的信號分類，包括預處理、特征提取和選擇、分類三個過程。國內外的眾多學者針對EEG信號分類過程做出了諸多探索，其中EEG融合特征提取和選擇是一個研究熱點。在EEG信號的帶通功率和功率譜密度(power spectral density,PSD)特征融合基礎上，Izabela等采用主成分分析(principal component analysis,PCA)算法[2]和遺傳算法(genetic algorithm,GA)[3]進行特征選擇；Coelho等[4]引入了一種新的人工免疫網絡算法cob-aiNet實現特征選擇；Bhattacharyya等[5]則針對特征選擇過程改進了微分進化算法，提出LA-DE算法；Atyabi等[6]應用遺傳算法和粒子群算法進行通道和特征選擇，均取得了較好的實驗結果。在EEG信號的共空間模式(common spatial patterns,CSP)特征基礎上，Sannelli等[7]提出了迭代通道選擇算法；Lin等[8]提出了結合Rayleigh系數的遺傳算法RC-GA算法用于通道選擇。Arvaneh等[9]引入稀疏共空間模式方法提高了分類準確率。Yuan等[10]基于有重疊的時間分段提取信號的F-score特征，有效降低了通道數量。German等[11]提出了同時采用PSD、Hjorth系數、自回歸(autoregressive,AR)模型參數和連續小波變換(continuous wavelet transform,CWT)4類特征提取方法的策略，引入最小角回歸(least angle regression,LARS)算法進行特征選擇。Pedro等[12]在提取帶通功率、Hjorth系數、自適應自回歸模型參數3種特征的基礎上采用局部Fisher判別分析(local Fisher discriminant analysis，LFDA)算法進行特征降維，通過特征變換(feature transformation,FT)提高了分類性能。文獻[11]和[12]均采用了提取多種EEG信號特征再降維處理的策略，同時提取多種EEG信號特征組成高維異構融合特征集合為特征選擇提供了更全面的候選特征，在此基礎上采用適當的特征選擇方法進行降維將有效提高分類性能。
聯合應用多種特征提取方法建立高維融合特征與特征降維的難點在于，全面的候選特征集合導致融合特征維數升高，而高維融合特征對特征降維方法提出了更高的要求。本文針對EEG信號的融合特征提取、最優特征子集選擇和分類提出了一種新方法。首先，聯合應用多種特征提取方法，建立高維融合特征；之后，采用封裝方式進行特征選擇：對訓練數據采用彈性網罰logistic回歸擬合模型，通過坐標下降法估計模型參數，運用10倍交叉驗證選擇出最優特征子集；最后采用已訓練的最優模型對測試樣本進行分類。
本文采用的特征提取方法是：時域統計、PSD、CSP和AR模型參數，這是目前應用最為廣泛的EEG信號特征提取方法。聯合應用這些方法，將提取的各種特征進行融合，建立高維特征集合是一種新的有益嘗試，這給后續的特征選擇提供了更全面的候選特征集合。對于特征選擇方法，本文對比了過濾式和封裝式兩種策略，包括彈性網罰邏輯回歸方法和Fisher準則下的前向特征選擇(forward feature selection,FFS)方法。訓練模型時采用了10倍交叉驗證方法優選特征子集。測試過程中對比了樸素貝葉斯分類器(naive Bayes classifier,NBC)、支持向量機(support vector machine,SVM)和邏輯回歸分類器。
1 特征提取
特征提取是EEG信號分類的關鍵技術之一，在結合神經生理學先驗知識的基礎上，本文聯合應用時域分析、頻域分析、空間分析和時空分析4類特征提取方法，建立高維融合特征。
時域統計量是一種常用的EEG信號特征，其提取方法簡單有效。本文實驗中對EEG信號的每個通道提取4個統計特征：原信號均值和標準差，以及原信號一階差分絕對值的均值和二階差分絕對值的均值。
不同運動想象任務激活的腦區不同，事件相關去同步(event related desynchronization,ERD)和事件相關同步(event related synchronization,ERS)與主動運動意識密切相關，出現在α節律(8~13 Hz)、β節律(13~30 Hz)和μ節律(8~12 Hz)的頻段，可以通過功率譜估計分析EEG信號節律的分布和變化[13]。實驗中對每個通道的信號提取5個頻段的帶通功率作為特征，分別是2~4、4~8、8~12、12~18和18~30 Hz頻段。
CSP方法是EEG信號特征提取的有效方法，對二分類情況的應用取得了顯著的成績[14]。該方法利用兩個協方差矩陣的同時對角化來設計一組空間濾波器，使兩類信號投影后能夠最大限度地被區分。實驗中針對不同通道數的一組EEG信號，均提取4維CSP特征。
AR模型是時間序列建模的有效工具，在BCI系統中得到了廣泛的應用[15]。實驗中對每個通道的信號建立6階AR模型，并將模型系數作為EEG信號特征。
2 特征選擇
聯合應用上述4類特征建立的融合特征必然是高維的，而訓練數據樣本卻不夠充分，這使分類器的訓練面臨高維小樣本的問題。同時，融合特征中不乏某些無關或冗余特征，直接輸入分類器進行分類不但增加了計算量，還會影響分類效果，這一結論也在實驗中得到了證實。因此，需要采用恰當的特征選擇方法，選出與分類最相關的特征子集。
實驗中采用了FFS方法和彈性網方法，FFS是一種簡單的自下而上的搜索方法，這里不再贅述，我們采用Fisher準則作為前向特征選擇的評價準則。
本文將彈性網罰邏輯回歸方法用于EEG信號特征選擇，通過極大化彈性網罰對數似然函數，得到模型的參數估計，即找到與分類相關的特征集合，完成選擇過程。其中，對目標函數的求解采用坐標下降法實現，該方法有效降低了計算量，使優化問題很快完成。
當響應變量是二元時，邏輯回歸模型經常使用。定義g表示響應變量，取值g={1,2}，預測變量xi=(xi1,xi2,…,xip)T，有N對觀測(xi,gi)，為簡化計算，假設xij已經標準化。
邏輯回歸模型通過預測變量的線性函數表示條件概率：
| $\begin{align} & \Pr (g=1x)=\frac{1}{1+{{e}^{-({{\beta }_{0}}+{{x}^{T}}\beta )}}}, \\ & \Pr (g=2x)=\frac{1}{1+{{e}^{+({{\beta }_{0}}+{{x}^{T}}\beta )}}} \\ \end{align}$ |
其中，β是p維的系數向量。用極大似然估計擬合上述模型，實驗數據集可以視為N次相互獨立的觀測，則有對數似然函數：
| $\begin{align} & L({{\beta }_{0}},\beta )=\frac{1}{N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,[{{y}_{i}}lnp({{x}_{i}})+(1-{{y}_{i}})ln(1-p({{x}_{i}}))] \\ & =\frac{1}{N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,[{{y}_{i}}({{\beta }_{0}}+x_{i}^{T}\beta )+ln(1-p({{x}_{i}}))] \\ \end{align}$ |
其中，。
在極大化對數似然函數的基礎上加入彈性網罰Pα(β)，則參數估計可以描述為：
| $({{\beta }_{0}},\beta )=\arg \underset{({{\beta }_{0}},\beta )\in {{R}^{P+1}}}{\mathop{\max }}\,\{L({{\beta }_{0}},\beta )-\lambda {{P}_{a}}(\beta )\}$ |
其中，。
可見彈性網罰是脊回歸罰和lasso罰的凸結合[16]，而lasso罰可以視為彈性網罰的一種特例，即當參數α=1時，彈性網罰變為lasso罰。
對式(3) 最大化問題的求解采用坐標下降法[17]。坐標下降法算法基本思想是將預測變量之間互不相關的多變量問題轉換為多個單變量的子問題，每次只優化一維變量，且優化系數可以在變量循環中更新，因此整個迭代過程將很快完成。
應用坐標下降法之前需要先做一項轉換。假設參數當前估計是{β0,β}，在當前估計值點做Taylor展開，可以得到式(2) 對數似然函數的二次逼近，即：
| ${{L}_{Q}}({{\beta }_{0}},\beta )=-\frac{1}{2N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}{{({{z}_{i}}-{{\beta }_{0}}-x_{i}^{T}\beta )}^{2}}+C{{({{\beta }_{0}},\beta )}^{2}}$ |
其中，視為響應；ωi=p(xi)(1-p(xi))為權值；，僅對βj進行部分優化時是常量；p(xi)是根據當前參數估計計算所得的值。
用式(4) 逼近形式替換式(3) 中的對數似然部分，則問題轉化為罰加權最小二乘形式：
| $({{\beta }_{0}},\beta )=arg\underset{({{\beta }_{0}},\beta )\in {{R}^{P+1}}}{\mathop{min}}\,\{-{{L}_{Q}}({{\beta }_{0}},\beta )+\lambda {{P}_{a}}(\beta )\}$ |
對式(5) 應用坐標下降法求解，記：
| $\begin{align} & -{{L}_{Q}}({{\beta }_{0}},\beta )+\lambda {{P}_{a}}(\beta )=\frac{1}{2N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}({{z}_{i}}-g{{({{x}_{i}})}^{(j)}} \\ & -{{x}_{ij}}{{\beta }_{j}}{{)}^{2}}+C{{({{\beta }_{0}},\beta )}^{2}}+\lambda {{P}_{a}}(\beta ) \\ \end{align}$ |
其中，，是除去xij一項的聯系函數。
假設參數當前估計是{β0,β}，每次只優化系數β的一維，其他維視為常數，則可以對系數β的第j維βj求導：
當βj＞0時，令導數等于0，可以求得βj的坐標更新形式：
| $\begin{matrix} {{\beta }_{j}}=\frac{\frac{1}{N}\underset{i-1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}{{x}_{ij}}({{z}_{i}}-g{{({{x}_{i}})}^{(j)}})-\lambda \alpha }{\frac{1}{N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}x_{ij}^{2}+\lambda (1-\alpha )} & ({{\beta }_{j}}>0) \\ \end{matrix}$ |
當βj＜0時，同樣可以得到類似的表達式，其他情況下βj=0，即有軟閾值算子描述的坐標更新形式：
| ${{\beta }_{j}}=\frac{S\left( \frac{1}{N}\underset{i-1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}{{x}_{ij}}({{z}_{i}}-g{{({{x}_{i}})}^{(j)}}),\lambda \alpha \right)}{\frac{1}{N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}x_{ij}^{2}+\lambda (1-\alpha )}$ |
總地來說，根據當前觀測值，建立彈性網罰邏輯回歸模型，采用坐標下降法求解系數β的過程是一系列循環的迭代求解過程，每個循環嵌套于上一個循環中，直到收斂。坐標下降法的步驟如下：
循環1：遞減的λ值；
循環2：使用當前參數{β0,β}更新LQ二次逼近；
循環3：通過坐標下降法求解式(5) 的罰加權最小二乘問題，得到βj，如式(8) 所示。
3 實驗過程與結果分析
3.1 數據說明
實驗中采用的EEG數據是國際BCI競賽Ⅳ的EEG數據，這組數據可以從網站下載，網址為http://www.bbci.de/competition/iv/#dataset1。
原始數據采用了59個電極通道，采樣頻率為1 000 Hz，每次實驗持續時間4 s，我們選用了其中11個電極的信號[18]進行分析，分別是FC3、FC4、Cz、C3、C4、C5、C6、T7、T8、CCP3和CCP4，并對相應信號重采樣，采樣頻率為100 Hz。預處理過程將原始數據的每次實驗部分對應的信號視為一個樣本，實驗任務視為該樣本的標簽，不考慮中間休息過程的信號，將原始信號整理為離散型樣本以便于分類訓練和測試。
3.2 實驗結果與分析
實驗過程首先對EEG信號進行特征提取，建立的融合特征維數為169維。對候選特征首先進行特征選擇，采用Fisher準則下的FFS方法從訓練集中選擇特征構成子集；之后將特征子集對應的訓練數據送入分類器進行100次10倍交叉驗證，得到該子集對應的平均錯誤率；針對不同特征數的子集，選擇平均錯誤率最小的作為最優子集；最后采用最優子集對應的訓練數據訓練分類器，將最優子集對應的測試數據送入分類器進行分類，得到測試錯誤率。以數據集A為例，不同特征數的子集對應的訓練錯誤率如圖 1所示。
圖1
數據集A融合特征FFS方法的訓練錯誤率
Figure1.
Train error of the fused feature on dataset A based on FFS method
對數據集A單獨應用時域統計量、帶通功率、CSP、AR模型系數4種特征提取方法以及聯合應用4種方法建立融合特征對應的實驗結果如表 2所示。實驗中特征選擇方法為FFS方法，預測過程分別采用支持向量機、邏輯回歸、貝葉斯3種分類器。不同分類器選出最優子集所需訓練時間如表 3所示。
表1
預處理后的實驗數據的說明
Table1.
Description of preprocessed data used in the experiment
表2
不同特征提取方法對應的測試錯誤率(%)
Table2.
Test error rate of different feature extraction methods (%)
表3
不同方法選出最優子集所需時間(s)
Table3.
Time required for different methods selecting the best subset (s)
可以看出，一方面，不經過特征選擇時，融合特征對EEG特征的提取雖然更全面，但特征維數高，特征存在冗余，導致分類效果并不好；經過特征選擇后，特征維數顯著減少，各種方法的正確率普遍提高，這證明了特征選擇過程的必要性。另一方面，經過特征選擇的融合特征的正確率雖然得到了提高，但效果并不明顯，且在訓練過程中，時間開銷隨著特征數的增加顯著變大，對100個子集逐個計算訓練正確率的過程耗時很長，對于高維數據不能滿足信號處理的實際需要。這說明FFS方法并不適合該問題，我們需要一個能夠自動優化子集特征數的快速特征選擇方法。
采用坐標下降法求解的彈性網罰邏輯回歸方法收斂速度快，且適用于高維數據的特征選擇問題，本文將彈性網方法引入BCI信號特征選擇和分類問題中，采用彈性網方法對高維異構融合特征進行特征選擇降維。
首先將彈性網罰邏輯回歸分類器錯誤率作為特征選擇的評價準則，給定100個不同的λ候選值和6個α候選值，即600種特征選擇的可能(對應的子集包含特征數可能相同也可能不同)，將訓練數據送入彈性網罰邏輯回歸分類器，進行10倍交叉驗證，得到不同特征子集對應的平均錯誤率，如圖 2所示。然后，選擇最小平均錯誤率對應的特征子集作為最優特征子集。最后，采用最優特征子集對應的訓練數據訓練預測分類器，將最優子集對應的測試數據送入預測分類器進行分類，得到測試錯誤率。圖 2中豎直虛線表示訓練錯誤率最低的子集對應的正則化參數和自由度取值。
圖2
數據集A融合特征彈性網方法的訓練錯誤率
Figure2.
Train error of the fused feature on dataset A based on elastic net method
特征選擇和預測過程可以采用過濾式方法或封裝式方法。過濾式彈性網方法的預測過程采用的分類器與特征選擇過程的分類器不相同，特征選擇過程中采用彈性網方法給出特征子集，預測過程的分類器則分別采用支持向量機、邏輯回歸、貝葉斯3種分類器；而封裝式彈性網方法的預測過程采用的分類器仍是彈性網分類器。但無論過濾式彈性網方法還是封裝式彈性網方法，特征子集的選擇過程是相同的。上述不同方法的實驗結果如表 4所示，彈性網方法選出最優子集所需時間及最優子集特征數如表 5所示。
表4
不同特征選擇和分類方法結果比較
Table4.
Comparison of feature selection and classification among different methods
表5
彈性網方法訓練時間和特征數
Table5.
Training time and number of feature of the elastic net method
可以看出，一方面，對于過濾式特征選擇，彈性網方法與FFS方法相比表現基本相當，但訓練時間短，能夠選出特征數較少的子集，即與分類更相關的子集，說明彈性網方法更適合高維EEG融合特征的最優子集選擇問題；另一方面，封裝式彈性網方法取得了更低的分類錯誤率，在預測準確性方面優于過濾式FFS特征選擇方法和過濾式彈性網方法。總地來說，封裝式彈性網方法對于EEG信號不僅分類效果好、速度快，還能給出相對簡單的回歸模型。
統計FFS方法結合支持向量機、邏輯回歸、貝葉斯3種分類器的特征選擇結果和彈性網方法的特征選擇結果包含的特征數，如圖 3所示。可以看出FFS結合貝葉斯分類器給出的特征子集包含特征數最少，其他方法表現相近，其中彈性網方法表現較為穩定。
圖3
不同方法的最優子集包含特征數
Figure3.
Feature number in the best subsets obtained using different methods
對于彈性網方法給出的最優子集進行分析，每一類特征在最優子集中所占的比例如圖 4所示。
圖4
所選特征在最優子集中所占比例
Figure4.
Proportion in selected features in each subject
融合特征全面地提取了EEG信號特征，特征選擇得到的最優子集中包含了全部4類特征，說明融合特征結合彈性網方法能夠將各種特征有機地結合起來，使各種特征提取方法聯合發揮作用，從而得到更有辨識能力的特征子集。
4 結論
信號分類是BCI系統的核心步驟，有效的特征提取和選擇方法是提高EEG信號識別正確率的關鍵。本文針對EEG信號分類，提出了一種封裝式彈性網新方法。首先，對預處理后的數據聯合應用多種特征提取方法，建立高維融合特征；其次，采用封裝式彈性網方法進行特征選擇：對訓練數據采用彈性網罰邏輯回歸模型擬合，通過坐標下降法求得模型的參數估計，運用10倍交叉驗證選出最優特征子集；最后對測試樣本采用已訓練的模型進行分類。其中，特征提取采用的方法是：時域統計特征、功率譜估計、CSP和AR模型系數。融合多種特征建立特征集合以供選擇是EEG信號分類方法的研究熱點。實驗證明，本文提出的封裝式彈性網罰邏輯回歸方法與其他方法相比，能夠更全面地提取EEG信號特征，選出與分類更相關的子集，提高了信號識別正確率，且算法更穩定、時間成本更低。對國際BCI競賽Ⅳ的第一組兩類運動想象數據的平均測試正確率達到了81.78%。下一步研究的問題是將封裝式彈性網特征選擇方法推廣到其它高維融合特征集合的子集選擇問題。
引言
腦機接口(brain-computer interface,BCI)是大腦和計算機或電子設備之間傳遞信息的通訊系統，這個系統不依賴外圍神經和肌肉組成的大腦正常通路[1]。這種全新的人機交互方式能夠將大腦產生的信息轉換成對外部設備的控制命令，從而代替肢體和語言實現與外界的交流。這項技術不僅能夠為思維正常但有運動障礙的人提供幫助，還能作為特殊環境下正常人對外部設備的控制方式，甚至能夠提供一種全新娛樂體驗。
BCI系統的核心問題之一是腦電圖(electroencephalogram,EEG)的信號分類，包括預處理、特征提取和選擇、分類三個過程。國內外的眾多學者針對EEG信號分類過程做出了諸多探索，其中EEG融合特征提取和選擇是一個研究熱點。在EEG信號的帶通功率和功率譜密度(power spectral density,PSD)特征融合基礎上，Izabela等采用主成分分析(principal component analysis,PCA)算法[2]和遺傳算法(genetic algorithm,GA)[3]進行特征選擇；Coelho等[4]引入了一種新的人工免疫網絡算法cob-aiNet實現特征選擇；Bhattacharyya等[5]則針對特征選擇過程改進了微分進化算法，提出LA-DE算法；Atyabi等[6]應用遺傳算法和粒子群算法進行通道和特征選擇，均取得了較好的實驗結果。在EEG信號的共空間模式(common spatial patterns,CSP)特征基礎上，Sannelli等[7]提出了迭代通道選擇算法；Lin等[8]提出了結合Rayleigh系數的遺傳算法RC-GA算法用于通道選擇。Arvaneh等[9]引入稀疏共空間模式方法提高了分類準確率。Yuan等[10]基于有重疊的時間分段提取信號的F-score特征，有效降低了通道數量。German等[11]提出了同時采用PSD、Hjorth系數、自回歸(autoregressive,AR)模型參數和連續小波變換(continuous wavelet transform,CWT)4類特征提取方法的策略，引入最小角回歸(least angle regression,LARS)算法進行特征選擇。Pedro等[12]在提取帶通功率、Hjorth系數、自適應自回歸模型參數3種特征的基礎上采用局部Fisher判別分析(local Fisher discriminant analysis，LFDA)算法進行特征降維，通過特征變換(feature transformation,FT)提高了分類性能。文獻[11]和[12]均采用了提取多種EEG信號特征再降維處理的策略，同時提取多種EEG信號特征組成高維異構融合特征集合為特征選擇提供了更全面的候選特征，在此基礎上采用適當的特征選擇方法進行降維將有效提高分類性能。
聯合應用多種特征提取方法建立高維融合特征與特征降維的難點在于，全面的候選特征集合導致融合特征維數升高，而高維融合特征對特征降維方法提出了更高的要求。本文針對EEG信號的融合特征提取、最優特征子集選擇和分類提出了一種新方法。首先，聯合應用多種特征提取方法，建立高維融合特征；之后，采用封裝方式進行特征選擇：對訓練數據采用彈性網罰logistic回歸擬合模型，通過坐標下降法估計模型參數，運用10倍交叉驗證選擇出最優特征子集；最后采用已訓練的最優模型對測試樣本進行分類。
本文采用的特征提取方法是：時域統計、PSD、CSP和AR模型參數，這是目前應用最為廣泛的EEG信號特征提取方法。聯合應用這些方法，將提取的各種特征進行融合，建立高維特征集合是一種新的有益嘗試，這給后續的特征選擇提供了更全面的候選特征集合。對于特征選擇方法，本文對比了過濾式和封裝式兩種策略，包括彈性網罰邏輯回歸方法和Fisher準則下的前向特征選擇(forward feature selection,FFS)方法。訓練模型時采用了10倍交叉驗證方法優選特征子集。測試過程中對比了樸素貝葉斯分類器(naive Bayes classifier,NBC)、支持向量機(support vector machine,SVM)和邏輯回歸分類器。
1 特征提取
特征提取是EEG信號分類的關鍵技術之一，在結合神經生理學先驗知識的基礎上，本文聯合應用時域分析、頻域分析、空間分析和時空分析4類特征提取方法，建立高維融合特征。
時域統計量是一種常用的EEG信號特征，其提取方法簡單有效。本文實驗中對EEG信號的每個通道提取4個統計特征：原信號均值和標準差，以及原信號一階差分絕對值的均值和二階差分絕對值的均值。
不同運動想象任務激活的腦區不同，事件相關去同步(event related desynchronization,ERD)和事件相關同步(event related synchronization,ERS)與主動運動意識密切相關，出現在α節律(8~13 Hz)、β節律(13~30 Hz)和μ節律(8~12 Hz)的頻段，可以通過功率譜估計分析EEG信號節律的分布和變化[13]。實驗中對每個通道的信號提取5個頻段的帶通功率作為特征，分別是2~4、4~8、8~12、12~18和18~30 Hz頻段。
CSP方法是EEG信號特征提取的有效方法，對二分類情況的應用取得了顯著的成績[14]。該方法利用兩個協方差矩陣的同時對角化來設計一組空間濾波器，使兩類信號投影后能夠最大限度地被區分。實驗中針對不同通道數的一組EEG信號，均提取4維CSP特征。
AR模型是時間序列建模的有效工具，在BCI系統中得到了廣泛的應用[15]。實驗中對每個通道的信號建立6階AR模型，并將模型系數作為EEG信號特征。
2 特征選擇
聯合應用上述4類特征建立的融合特征必然是高維的，而訓練數據樣本卻不夠充分，這使分類器的訓練面臨高維小樣本的問題。同時，融合特征中不乏某些無關或冗余特征，直接輸入分類器進行分類不但增加了計算量，還會影響分類效果，這一結論也在實驗中得到了證實。因此，需要采用恰當的特征選擇方法，選出與分類最相關的特征子集。
實驗中采用了FFS方法和彈性網方法，FFS是一種簡單的自下而上的搜索方法，這里不再贅述，我們采用Fisher準則作為前向特征選擇的評價準則。
本文將彈性網罰邏輯回歸方法用于EEG信號特征選擇，通過極大化彈性網罰對數似然函數，得到模型的參數估計，即找到與分類相關的特征集合，完成選擇過程。其中，對目標函數的求解采用坐標下降法實現，該方法有效降低了計算量，使優化問題很快完成。
當響應變量是二元時，邏輯回歸模型經常使用。定義g表示響應變量，取值g={1,2}，預測變量xi=(xi1,xi2,…,xip)T，有N對觀測(xi,gi)，為簡化計算，假設xij已經標準化。
邏輯回歸模型通過預測變量的線性函數表示條件概率：
| $\begin{align} & \Pr (g=1x)=\frac{1}{1+{{e}^{-({{\beta }_{0}}+{{x}^{T}}\beta )}}}, \\ & \Pr (g=2x)=\frac{1}{1+{{e}^{+({{\beta }_{0}}+{{x}^{T}}\beta )}}} \\ \end{align}$ |
其中，β是p維的系數向量。用極大似然估計擬合上述模型，實驗數據集可以視為N次相互獨立的觀測，則有對數似然函數：
| $\begin{align} & L({{\beta }_{0}},\beta )=\frac{1}{N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,[{{y}_{i}}lnp({{x}_{i}})+(1-{{y}_{i}})ln(1-p({{x}_{i}}))] \\ & =\frac{1}{N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,[{{y}_{i}}({{\beta }_{0}}+x_{i}^{T}\beta )+ln(1-p({{x}_{i}}))] \\ \end{align}$ |
其中，。
在極大化對數似然函數的基礎上加入彈性網罰Pα(β)，則參數估計可以描述為：
| $({{\beta }_{0}},\beta )=\arg \underset{({{\beta }_{0}},\beta )\in {{R}^{P+1}}}{\mathop{\max }}\,\{L({{\beta }_{0}},\beta )-\lambda {{P}_{a}}(\beta )\}$ |
其中，。
可見彈性網罰是脊回歸罰和lasso罰的凸結合[16]，而lasso罰可以視為彈性網罰的一種特例，即當參數α=1時，彈性網罰變為lasso罰。
對式(3) 最大化問題的求解采用坐標下降法[17]。坐標下降法算法基本思想是將預測變量之間互不相關的多變量問題轉換為多個單變量的子問題，每次只優化一維變量，且優化系數可以在變量循環中更新，因此整個迭代過程將很快完成。
應用坐標下降法之前需要先做一項轉換。假設參數當前估計是{β0,β}，在當前估計值點做Taylor展開，可以得到式(2) 對數似然函數的二次逼近，即：
| ${{L}_{Q}}({{\beta }_{0}},\beta )=-\frac{1}{2N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}{{({{z}_{i}}-{{\beta }_{0}}-x_{i}^{T}\beta )}^{2}}+C{{({{\beta }_{0}},\beta )}^{2}}$ |
其中，視為響應；ωi=p(xi)(1-p(xi))為權值；，僅對βj進行部分優化時是常量；p(xi)是根據當前參數估計計算所得的值。
用式(4) 逼近形式替換式(3) 中的對數似然部分，則問題轉化為罰加權最小二乘形式：
| $({{\beta }_{0}},\beta )=arg\underset{({{\beta }_{0}},\beta )\in {{R}^{P+1}}}{\mathop{min}}\,\{-{{L}_{Q}}({{\beta }_{0}},\beta )+\lambda {{P}_{a}}(\beta )\}$ |
對式(5) 應用坐標下降法求解，記：
| $\begin{align} & -{{L}_{Q}}({{\beta }_{0}},\beta )+\lambda {{P}_{a}}(\beta )=\frac{1}{2N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}({{z}_{i}}-g{{({{x}_{i}})}^{(j)}} \\ & -{{x}_{ij}}{{\beta }_{j}}{{)}^{2}}+C{{({{\beta }_{0}},\beta )}^{2}}+\lambda {{P}_{a}}(\beta ) \\ \end{align}$ |
其中，，是除去xij一項的聯系函數。
假設參數當前估計是{β0,β}，每次只優化系數β的一維，其他維視為常數，則可以對系數β的第j維βj求導：
當βj＞0時，令導數等于0，可以求得βj的坐標更新形式：
| $\begin{matrix} {{\beta }_{j}}=\frac{\frac{1}{N}\underset{i-1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}{{x}_{ij}}({{z}_{i}}-g{{({{x}_{i}})}^{(j)}})-\lambda \alpha }{\frac{1}{N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}x_{ij}^{2}+\lambda (1-\alpha )} & ({{\beta }_{j}}>0) \\ \end{matrix}$ |
當βj＜0時，同樣可以得到類似的表達式，其他情況下βj=0，即有軟閾值算子描述的坐標更新形式：
| ${{\beta }_{j}}=\frac{S\left( \frac{1}{N}\underset{i-1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}{{x}_{ij}}({{z}_{i}}-g{{({{x}_{i}})}^{(j)}}),\lambda \alpha \right)}{\frac{1}{N}\underset{i=1}{\overset{N}{\mathop{\sum }}}\,{{\omega }_{i}}x_{ij}^{2}+\lambda (1-\alpha )}$ |
總地來說，根據當前觀測值，建立彈性網罰邏輯回歸模型，采用坐標下降法求解系數β的過程是一系列循環的迭代求解過程，每個循環嵌套于上一個循環中，直到收斂。坐標下降法的步驟如下：
循環1：遞減的λ值；
循環2：使用當前參數{β0,β}更新LQ二次逼近；
循環3：通過坐標下降法求解式(5) 的罰加權最小二乘問題，得到βj，如式(8) 所示。
3 實驗過程與結果分析
3.1 數據說明
實驗中采用的EEG數據是國際BCI競賽Ⅳ的EEG數據，這組數據可以從網站下載，網址為http://www.bbci.de/competition/iv/#dataset1。
原始數據采用了59個電極通道，采樣頻率為1 000 Hz，每次實驗持續時間4 s，我們選用了其中11個電極的信號[18]進行分析，分別是FC3、FC4、Cz、C3、C4、C5、C6、T7、T8、CCP3和CCP4，并對相應信號重采樣，采樣頻率為100 Hz。預處理過程將原始數據的每次實驗部分對應的信號視為一個樣本，實驗任務視為該樣本的標簽，不考慮中間休息過程的信號，將原始信號整理為離散型樣本以便于分類訓練和測試。
3.2 實驗結果與分析
實驗過程首先對EEG信號進行特征提取，建立的融合特征維數為169維。對候選特征首先進行特征選擇，采用Fisher準則下的FFS方法從訓練集中選擇特征構成子集；之后將特征子集對應的訓練數據送入分類器進行100次10倍交叉驗證，得到該子集對應的平均錯誤率；針對不同特征數的子集，選擇平均錯誤率最小的作為最優子集；最后采用最優子集對應的訓練數據訓練分類器，將最優子集對應的測試數據送入分類器進行分類，得到測試錯誤率。以數據集A為例，不同特征數的子集對應的訓練錯誤率如圖 1所示。
圖1
數據集A融合特征FFS方法的訓練錯誤率
Figure1.
Train error of the fused feature on dataset A based on FFS method
對數據集A單獨應用時域統計量、帶通功率、CSP、AR模型系數4種特征提取方法以及聯合應用4種方法建立融合特征對應的實驗結果如表 2所示。實驗中特征選擇方法為FFS方法，預測過程分別采用支持向量機、邏輯回歸、貝葉斯3種分類器。不同分類器選出最優子集所需訓練時間如表 3所示。
表1
預處理后的實驗數據的說明
Table1.
Description of preprocessed data used in the experiment
表2
不同特征提取方法對應的測試錯誤率(%)
Table2.
Test error rate of different feature extraction methods (%)
表3
不同方法選出最優子集所需時間(s)
Table3.
Time required for different methods selecting the best subset (s)
可以看出，一方面，不經過特征選擇時，融合特征對EEG特征的提取雖然更全面，但特征維數高，特征存在冗余，導致分類效果并不好；經過特征選擇后，特征維數顯著減少，各種方法的正確率普遍提高，這證明了特征選擇過程的必要性。另一方面，經過特征選擇的融合特征的正確率雖然得到了提高，但效果并不明顯，且在訓練過程中，時間開銷隨著特征數的增加顯著變大，對100個子集逐個計算訓練正確率的過程耗時很長，對于高維數據不能滿足信號處理的實際需要。這說明FFS方法并不適合該問題，我們需要一個能夠自動優化子集特征數的快速特征選擇方法。
采用坐標下降法求解的彈性網罰邏輯回歸方法收斂速度快，且適用于高維數據的特征選擇問題，本文將彈性網方法引入BCI信號特征選擇和分類問題中，采用彈性網方法對高維異構融合特征進行特征選擇降維。
首先將彈性網罰邏輯回歸分類器錯誤率作為特征選擇的評價準則，給定100個不同的λ候選值和6個α候選值，即600種特征選擇的可能(對應的子集包含特征數可能相同也可能不同)，將訓練數據送入彈性網罰邏輯回歸分類器，進行10倍交叉驗證，得到不同特征子集對應的平均錯誤率，如圖 2所示。然后，選擇最小平均錯誤率對應的特征子集作為最優特征子集。最后，采用最優特征子集對應的訓練數據訓練預測分類器，將最優子集對應的測試數據送入預測分類器進行分類，得到測試錯誤率。圖 2中豎直虛線表示訓練錯誤率最低的子集對應的正則化參數和自由度取值。
圖2
數據集A融合特征彈性網方法的訓練錯誤率
Figure2.
Train error of the fused feature on dataset A based on elastic net method
特征選擇和預測過程可以采用過濾式方法或封裝式方法。過濾式彈性網方法的預測過程采用的分類器與特征選擇過程的分類器不相同，特征選擇過程中采用彈性網方法給出特征子集，預測過程的分類器則分別采用支持向量機、邏輯回歸、貝葉斯3種分類器；而封裝式彈性網方法的預測過程采用的分類器仍是彈性網分類器。但無論過濾式彈性網方法還是封裝式彈性網方法，特征子集的選擇過程是相同的。上述不同方法的實驗結果如表 4所示，彈性網方法選出最優子集所需時間及最優子集特征數如表 5所示。
表4
不同特征選擇和分類方法結果比較
Table4.
Comparison of feature selection and classification among different methods
表5
彈性網方法訓練時間和特征數
Table5.
Training time and number of feature of the elastic net method
可以看出，一方面，對于過濾式特征選擇，彈性網方法與FFS方法相比表現基本相當，但訓練時間短，能夠選出特征數較少的子集，即與分類更相關的子集，說明彈性網方法更適合高維EEG融合特征的最優子集選擇問題；另一方面，封裝式彈性網方法取得了更低的分類錯誤率，在預測準確性方面優于過濾式FFS特征選擇方法和過濾式彈性網方法。總地來說，封裝式彈性網方法對于EEG信號不僅分類效果好、速度快，還能給出相對簡單的回歸模型。
統計FFS方法結合支持向量機、邏輯回歸、貝葉斯3種分類器的特征選擇結果和彈性網方法的特征選擇結果包含的特征數，如圖 3所示。可以看出FFS結合貝葉斯分類器給出的特征子集包含特征數最少，其他方法表現相近，其中彈性網方法表現較為穩定。
圖3
不同方法的最優子集包含特征數
Figure3.
Feature number in the best subsets obtained using different methods
對于彈性網方法給出的最優子集進行分析，每一類特征在最優子集中所占的比例如圖 4所示。
圖4
所選特征在最優子集中所占比例
Figure4.
Proportion in selected features in each subject
融合特征全面地提取了EEG信號特征，特征選擇得到的最優子集中包含了全部4類特征，說明融合特征結合彈性網方法能夠將各種特征有機地結合起來，使各種特征提取方法聯合發揮作用，從而得到更有辨識能力的特征子集。
4 結論
信號分類是BCI系統的核心步驟，有效的特征提取和選擇方法是提高EEG信號識別正確率的關鍵。本文針對EEG信號分類，提出了一種封裝式彈性網新方法。首先，對預處理后的數據聯合應用多種特征提取方法，建立高維融合特征；其次，采用封裝式彈性網方法進行特征選擇：對訓練數據采用彈性網罰邏輯回歸模型擬合，通過坐標下降法求得模型的參數估計，運用10倍交叉驗證選出最優特征子集；最后對測試樣本采用已訓練的模型進行分類。其中，特征提取采用的方法是：時域統計特征、功率譜估計、CSP和AR模型系數。融合多種特征建立特征集合以供選擇是EEG信號分類方法的研究熱點。實驗證明，本文提出的封裝式彈性網罰邏輯回歸方法與其他方法相比，能夠更全面地提取EEG信號特征，選出與分類更相關的子集，提高了信號識別正確率，且算法更穩定、時間成本更低。對國際BCI競賽Ⅳ的第一組兩類運動想象數據的平均測試正確率達到了81.78%。下一步研究的問題是將封裝式彈性網特征選擇方法推廣到其它高維融合特征集合的子集選擇問題。
