CN102631194B - 一种用于心电特征选择的禁忌搜索方法 - Google Patents

Abstract

本发明公开了一种用于心电特征选择的禁忌搜索方法，包括如下步骤：采用序列后向算法，形成一个N-1行，N列的二维表L，其中N为选择的特征总维数，每一列代表一个特征，每一行称为一个空间，其中第n个空间选择有n个特征，1≤n≤N-1；表内每个元素的值用“0”或“1”表示，“0”代表在进行特征选择时该元素没有被选中，“1”代表该元素被选中；对每一空间内的选中的特征采用禁忌搜索算法进行求解，得到每个空间的解组成的表S；选择各个空间中适应度函数最大的作为最终特征选择结果。本发明采用在禁忌搜索算法中采用划分空间的方法结合Fisher分类器完成选择过程，克服了基本禁忌搜索算法容易陷入局部最优的缺点。

Description

一种用于心电特征选择的禁忌搜索方法

技术领域

本发明涉及一种生理信号情感识别方法，尤其涉及一种心电特征选择方法。

背景技术

情感识别是情感计算领域中的一个重要组成部分，如果要让计算机能从面部、姿势、语音和生理信号中识别出人类的情感，首先必须让计算机能识别情感，所以情感识别就显得非常重要。情感识别主要包括肢体情感识别、面部表情情感识别和语音情感识别和生理信号情感识别等方面。在“心率变异性在心身疾病和情绪障碍研究中的应用(心理学进展，2006，14(2)：261-265)”一文中，提出面部表情、语言姿势等流露出来的情感状态主要受神经系统控制，容易受到主观意识的控制。而生理信号是伴随着人的情感变化由人体内部器官产生的一种生物电信号，不易受主观意识的控制，所以更能客观真实的反映人的情感状态，更有鲁棒性和客观性，但基于生理信号情感识别的研究最困难的。

心电(ECG)信号是重要的生理信号之一，在“Affective PatternClassification(Perceptual Computing Section Technical Report，1998，473)”一文中，美国MIT实验室的Picard教授采用四种生理信号(EMG、ECG、GSR、RSP)率先证明从生理信号中提取特征模式来进行情感识别是有效的，并且研究表明，ECG信号能够体现情感状态的变化，从中提取的特征能够反映情感状态的差异。由于心电、心率信号的研究技术相当成熟，而且在人们在经历情感变化时，心电信号会有不同程度的变化趋势，因此，以心电信号作为研究对象进行情感识别的学者举不胜举。在“Affective Wearable(In Proceedings of the FirstInternational Symposium on Wearable Computers，Oct，1997，Cambridge，MA)”一文中，美国的麻省理工学院媒体实验室曾经从心电信号P-QRS-T波各波的间隔、幅度等计算六种统计特征(均值、中值、方差、最大值、最小值和范围)采用SFFS、Fisher投影特征选择算法和DFA、QDF分类器对8种情感进行分类，取得了较好的识别结果。在“Emotion Recognition System using Short-term Monitoring ofPhysiological Signals(Medical Biology Engine Computer，2004，42：419-427)”一文中，韩国的K.H.Kim采集175个被试的情感心电信号通过一定的设计方法，直接从心电信号提取出心率信号以及心率变异率信号(HRV)，从中提取特征后直接用支持向量机(SVM)对悲伤、愤怒、压力和惊奇四种情感进行分类；德国奥森堡大学的多媒体与信号处理实验室也是主要从心电信号提取出特征，着重比较了用SBS、SFS、ANOA这些特征选择方法和KNN、LDF、MLP不同分类器相结合的情感识别效果，并研究发现：较低的心率变异率(HRV)表明是放松的状态，而增强的HRV表明可能是精神紧张和受到挫折的状态。

在特征选择中，研究者大多采用传统的特征选择方法，有得甚至直接将提取的特征用于情感分类；而分类器的选择也是根据经验选择分类效果好的分类器，所以各种各样的分类器也就被应用在情感分类中。对特征的提取也多采用提取统计特征的方法。

在特征选择中，采用传统的禁忌搜索算法，常常会引发维数灾难，导致特征选择不能覆盖在所有可能的解，使得最终特征选择出的结果不为最优解。

发明内容

本发明的目的是提供一种有效避免维数灾难的用于心电特征选择中禁忌搜索方法。

本发明为实现上述目的，采用如下技术方案：一种用于心电特征选择的禁忌搜索方法，其特征在于：所述方法包括如下步骤：

采用序列后向算法，形成一个N-1行，N列的二维表L，其中N为选择的特征总维数，每一列代表一个特征，每一行称为一个空间，其中第n个空间选择有n个特征，1≤n≤N-1；表内每个元素的值用“0”或“1”表示，“0”代表在进行特征选择时该元素没有被选中，“1”代表该元素被选中；

对每一空间内的选中的特征采用禁忌搜索算法进行求解，得到每个空间的解组成的表S；

选择各个空间中适应度函数最大的作为最终特征选择结果。

其中每个空间中禁忌搜索算法的步骤为：

S1：初始设禁忌表T＝Ф，设置禁忌长度，设置最大迭代步数；将此空间内通过序列后向算法得到的值作为初始解，并把它作为该空间的暂定全局最优解Bestsofar和迭代搜索的起点，即当前局部最优解cand，计算初始解的适应度函数

值；其中ratA和ratB分别代表目标情感的正确识别率和非目标情感的正确识别率；

S2：判断是否满足最大迭代步数的停止准则，如满足则终止该空间的计算，将把Bestsofar放入表S中，结束此空间计算；

如不满足，则将当前局部最优解cand为下次的迭代起点；

S3：生成N个与该空间具有相同特征选择个数的候选解集；

S4：寻优：计算出每个候选解的适应度函数f的值，从候选解集中选出适应度函数值最大的解，

将该解与禁忌表中的解比较；如果该解不在禁忌表中，将其中较大的值作为当前局部最优解和Bestsofar的值；则转S5；

如果该解在禁忌表中，且不满足特赦准则，将适应度函数值第二大的解作为局部最优解；

如果该解在禁忌表中，且满足特赦准则，将禁忌表中的该解提前释放，并作为局部最优解，并将该解的适应度函数值与Bestsofar的适应度函数值比较，如大于，则用该解作为Bestsofar；其中特赦准则指的是，禁忌表中的某个值在迭代过程中作为局部最优解出现规定次数；

S5：更新禁忌表：将当前局部最优解写入禁忌表，转S2。

所述S3步骤中生成候选解集的方法为：分别将cand的第i位cand(i)值变为1-cand(i)，从第i位的下一位开始查找，遇到值等于1-cand(i)的p位，将其值变为1-cand(p)后结束；若P＝N，则P从1开始继续查，遇到值等于1-cand(i)的p位，将其值变为1-cand(p)后结束；

其中1≤i≤N，1≤P≤N，N为特征总维数。

本发明采用在禁忌搜索算法中采用划分空间的方法结合Fisher分类器完成选择过程，克服了基本禁忌搜索算法容易陷入局部最优的缺点。实验仿真结果表明，将相关性分析和改进的智能优化算法用于生理信号的情感识别领域是完全可行的。

说明书附图

图1为Fisher分类器中中间数据分布图；

图2为Fisher分类器中最终测试数据分布图。

具体实施方式

下面结合附图和具体实施例对本发明作进一步的阐述。

一种基于心电信号的情感识别特征选择方法，包括下列步骤：

1.情感心电信号数据采集流程；

实验设备采用美国Biopac公司提供的多导生理记录仪MP150。实验对象为年龄17-20岁本科一年级的391个在校学生，所有预约被试自愿参加本实验，身体健康，无心脏疾病和精神疾病病史。实验素材是从大量电影中挑选出的分别对高兴、惊奇、厌恶、悲伤、愤怒、恐惧六种情感具有很好诱发效果的片段。实验前，主试要让被试填写个人信息、自愿测试的同意书和一份测试被试是否有述情障碍的问卷调查。实验过程中，每个情感诱发片段后会让被试填写问卷调查来测试其是否有该种情感，并用5个尺度(如1代表“很弱”，5代表“很强”)来对情感强度作出评价。

2.ECG信号去噪；

在完成数据的采集后能够得到原始数据，心电信号的采集中，有来自自身的干扰如肌电干扰，也有来自外界的干扰如工频干扰等，从而导致信号的信噪比下降。而且不同被试针对同一个素材的ECG信号不一致，存在个体差异性；针对上述问题对采集的ECG数据进行去噪、标准化和归一化。

(2-1)去噪

采用小波变换的方法对心电信号进行P-QRS-T波检测，主要是通过合适的小波函数对心电信号进行连续小波变换来准确定位R波，在R波定位后相继检测Q波、S波、P波和T波。

小波变换的含义是：将基本小波ψ(t)作位移τ，在不同尺度α下与信号f(t)作内积：

${\mathrm{WT}}_f(a,\mathrm{\τ})=\frac{1}{\sqrt{a}}{\∫}_{-\∞}^{+\∞}f\left(t\right){\mathrm{\ψ}}^{*}\left(\frac{t-\mathrm{\τ}}{a}\right)\mathrm{dt},\mathrm{\α}>0---\left(1\right)$

其在频域的等效表示为：

${\mathrm{WT}}_f(a,\mathrm{\τ})=\frac{\sqrt{\mathrm{\α}}}{2\mathrm{\π}}\∫F\left(\mathrm{\ω}\right){\mathrm{\ψ}}^{*}\left(\mathrm{\α\ω}\right)e^{\mathrm{j\ω\τ}}\mathrm{d\ω},\mathrm{\α}>0---\left(2\right)$

其中，F(ω)，ψ(ω)分别是f(t)，ψ(t)的傅立叶变换。小波变换相当于把ψ(t)当作镜头对信号f(t)放大或缩小。τ的作用相当于让镜头对目标平行移动，α的作用相当于镜头向目标推进或远离。

由于强噪声干扰会引起虚假R峰，从而导致异常R-R间隔出现。文中通过异常R-R间隔作为强噪声干扰的标志，避免相应信号段的特征提取。Q、S、P、T波根据R波位置从原始心电信号中定位。一个心电周期中，Q、S波对应于该周期R波位置左，右波谷位置。P、T波需要在相邻QRS波群之间经过3次几何平均平滑降噪处理后定位。消除高频噪声后，分别检测R-R间隔17％-50％和75％-83.3％间的最大值处即为T波和P波位置。

在相邻QRS波群直接采用几何平均进行高频降噪处理，其几何平均公式(3)如下，x(n)表示原始心电信号。

$x\left(n\right)=\frac{x(n-1)+x(n+1)}{2}---\left(3\right)$

下面举例来说明去噪过程：

首先对原始心电信号进行采样，以恐惧信号为例，对该信号截取了80秒，合计16000个采样点的情感生理信号，n的最大值为16000。

对一个被试采集得到的心电信号进行抽样的第一个采样点值为790.15，第三个采样点值为792.16，然后降噪后第二个采样点值为(790.15+792.17)/2＝791.16，然后依次类推进行平滑，得到经过降噪后的恐惧的心电信号。

(2-2)数据的标准化：由于招募的被试者不同，其ECG信号存在着个体差异性，这种差异性体现在不同的人在同一个时间面对相同的环境，以及同一个人在不同时间面对不同环境下得出的。为了建立通用的基于ECG信号的情感识别系统，必须去除这种个体差异性后，得到的模型才能被推广，具体的操作如式4所示：

D_o＝D_emotion-D_calm    (4)

其中，D_emotion为降噪后的原始数据，D_calm为同一个被试者在由素材中的风景画和轻音乐环境下记录的ECG信号数据，得到的D₀即为经过标准化之后的数据。

(2-3)归一化处理

提取数据统计特征时，情感特征取值用数据基线的相应统计量进行归一化。

具体公式如下：

设其中一个特征数据为x＝(x₁，x₂，.....，x_m)，m表示被试人数，x为特征数据，建立映射f如下：

$f\left(k\right)=\frac{x_k-x_{\min }}{x_{\max }-x_{\min }}$

公式中的参数表示为x_max＝max(x)＝max(x₁，x₂，...，x_m)，x_min＝min(x)＝min(x₁，x₂，...，x_m)，x_k是原始特征数据。

下面以恐惧情感例举一个样本数据的具体处理过程如下：

提取恐惧的1号特征(‘ecg-mean’)的一个被试的值为785.9091，在1号特征的204个被试的数值值中最大值为1165.7，最小值为501.3782，然后将该被试的1号特征值归一化值＝(785.9091-501.3782)/(1165.7-501.3782)＝0.4283

3.ECG信号特征提取

在准确定位P-QRS-T波后，从中提取可能反映情感变化的111个统计特征作为原始特征集构成原始特征集合。其中包括：相邻P、Q、R、S和T波间隔，P-Q、Q-S、S-T间隔，P、Q、R、S和T波峰值，相邻R、T波幅度的一阶差分以及它们的一阶差分的绝对值，QRS波群、P波和T波各波能量的均值，中值，标准差，最小值，最大值，最大范围(最大与最小值之差)。采用在“基于LabVIEW的小波变换在心电信号分析中的应用研究(合肥：合肥工业大学，2008.10-20)”一硕士论文中，提到的Daubechies6(db6)小波对原始心电信号做离散小波变换8层分解。将第3和4层的细节系数进行重构，计算其QRS波群位置的重构信号的平均值和方差分别占整个重构信号的比值；同时，将第3、4和5层的细节系数进行重构，计算其QRS波群和T位置的重构信号的平均值和方差分别占整个重构信号的比值，并将其比值作为特征。

心电信号的主要特征提取公式：

1)均值

${\mathrm{\μ}}_x=\frac{1}{N}\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{x}_n---\left(5\right)$

2)均方差

${\mathrm{\σ}}_x={\left(\frac{1}{N-1}\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{(\mathrm{Xn}-\mathrm{\μx})}^2\right)}^{1/2}---\left(6\right)$

3)一阶差分绝对值的均值

${\mathrm{Diff}}_1\left(x\right)=\frac{1}{N-1}\underset{n=1}{\overset{N-1}{\mathrm{\Σ}}}|X_{n+1}-X_n|---\left(7\right)$

4)二阶差分绝对值的均值

${\mathrm{Diff}}_{2}f\left(x\right)=\frac{1}{N-2}\underset{n=1}{\overset{N-2}{\mathrm{\Σ}}}|X_{n+2}-X_n|---\left(8\right)$

5)最大最小值之差

range_i＝Max-Min    (9)

7)一阶差分的均方根

$\mathrm{rSD}=\sqrt{\frac{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(X_n-X_{n-1})}^2}{N}}---\left(10\right)$

4.相关性分析

在“Optimal discriminant analysis using the rank decomposition(Pattern Recognition，1992，25(1)：101-111)”一文中，Cheng、Zhang和Yang讨论了最佳鉴别矢量的维数与分类能力的关系，从理论上和实验结果分别证明出在原特征集合中增加或者删除相关特征并不改变最优特征子集的分类能力。

该方法用公式(11)计算出111个原始特征两两特征的相关系数，(r_xy表示线性相关系数；X，Y表示其中两个特征；x_i为X第i个样本值；y_i为Y第i个样本值；

为X均值；

为Y均值；n表示特征集的样本个数)

$r_{\mathrm{xy}}=|\left(\frac{1}{n}\underset{i}{\overset{n}{\mathrm{\Σ}}}(x_i-\stackrel{\‾}{x})(y_i-\stackrel{\‾}{y})\right)/\left(\sqrt{\frac{1}{n}\underset{i}{\overset{n}{\mathrm{\Σ}}}(x_i-\stackrel{\‾}{x})}\sqrt{\frac{1}{n}\underset{i}{\overset{n}{\mathrm{\Σ}}}(y_i-\stackrel{\‾}{y})}\right)|---\left(11\right)$

得到一个相关系数矩阵后，采用单独最优特征组合排序方法对111个特征对分类能力贡献大小进行排序，将r_xy＞δ(0≤δ≤1)的两个特征中删除一个排序相对靠后的特征，但当遇到两个特征中已有一个特征被删除时，要保留与被删除特征相关的哪个特征。

相关性分析步骤如下：

Step 1：删除原始特征集M_m×n(m是样本个数，n是特征维数为111)中畸形的样本，并且将特征矩阵归一化，得到特征矩阵M_l×n(l为删除样本后的样本个数)。

Step 2：采用在“特征选择算法研究综述(模式识别与人工智能，2007，4：V01.20 No.2)”一文中记录的后向算法计算111个原始特征在分别以Anger、Disgust、Fear、Grief、Joy和Surprise为目标情感状态下的适应度函数值，并结合Fisher分类器，逐个度量单个特征的分类能力，以此对特征矩阵M_l×n进行排序。

Step 3：计算出经过Step 1处理后的111个特征之间的线性相关系数矩阵。

Step 4：假定一个阀值δ(0≤δ≤1)，如果r_xy＞δ(0≤δ≤1)，参照Step 2的结果，保留最优的特征，剔除其余特征，将保留下来的特征构成一个降维后的特征矩阵M_l×p(p为经过删除后的特征维数)。设置不同的门限值δ可调节入选特征的数目，门限值越小，入选特征个数越少，特征间的相关性也越小。

这里举个4个特征的实例来说明该处理过程：

Step1：去除畸形样本和对样本进行标准化后得到4个样本数据如下：

表1采集的样本数据

Step2：对以上四个样本进行归一化，计算过程如下：

例如：1号特征的第一个样本数据进行归一化：

其它计算方法相同，得到归一化后的表如下：

表2归一化后的样本数据

Step3：将单个特征送入后向算法结合Fisher分类器中进行计算每个特征的分类识别率，然后采用识别率从大到小的方式对特征进行排序如下：

表3采用后向算法对特征贡献从大到小的排序

特征贡献排序

4

2

3

1

Step4：得到相关系数表为

表4特征之间的相关性系数表

特征号	1	2	3	4
					1	1.00	0.30	0.47	0.92
2	0.30	1.00	0.45	0.37
					3	0.47	0.45	1.00	0.09

4

0.92

0.37

0.09

1.00

Step5：这里取阀值δ(0≤δ≤1)为0.91，可看出1号特征和4号特征相关程度为0.92，大于0.91，所以要在1号和4号特征中删除一个特征，这时参考排序表，应该删除1号特征。

经过相关性分析后留下了2、3、4号特征。

5.基于改进的禁忌搜索算法进行心电信号的特征选择

禁忌搜索(Tabu Search，简称TS)是对人类思维过程本身的一种模拟，它通过对一些局部最优解的禁忌(也可以说是记忆)达到接纳一部分较差解，从而达到跳出局部搜索的目的。其最重要的思想是采用了禁忌技术，同时为了尽可能不错过产生最优解的“移动”，禁忌搜索还采用“藐视准则”策略，则可以无视其禁忌属性而仍采纳其为当前选择，以避免遗失优良状态。

该算法与传统的优化算法相比具有很好的爬山能力且区域集中搜索与全局分散搜索能较好平衡，为克服基本算法的不足，提出了一种改进的禁忌搜索算法应用于肌电信号的情感识别。该算法将搜索空间分解分为多个子空间，利用序列后向算法在每个子空间得到的解作为禁忌搜索的初始解，分类时，选用“图像模式识别(杨淑莹，清华大学出版社，2006，95-101)”一书中记录的线性判别函数Fisher分类器，适应度函数定义成分类正确识别率，并将需要分出来的一类情感看作目标情感，把其余五类情感看作非目标情感，如要对高兴与其它五类情感进行分类，则高兴为目标情感，其它五类情感为非目标情感，然后从每个初始解出发在对应的子空间开始搜索，得到每个子空间的最优特征子集，最后对这些解进行评价，找出整个空间的最优特征子集。

5.1适应度函数

评价得到的特征子集组合的优劣的依据是选用“图像模式识别(杨淑莹，清华大学出版社，2006，95-101)”一书中记录的线性判别函数Fisher分类器计算出的两类正确识别率，将适应度函数定义如下所示：

$f=\sqrt{({\mathrm{ratA}}^2+{\mathrm{ratB}}^2)}$

其中设ratA和ratB分别代表目标情感的正确识别率和非目标情感的正确识别率。下面举例说明ratA和ratB值在Fisher线性分类器中的计算过程。

Fisher分类的基本原理就是找到一个最合适的投影轴，使两类样本在该轴上的类间距尽可能大，样本类内间距尽可能小，这样可以达到最佳的分类效果，如何找到这个投影轴，就是要解决的首要问题。

Step 1：设置包含目标和非目标情感的样本数据集，用于求出最优的投影轴。设置数据如下：

目标情感包括4个样本和2个特征，数据如下：

$w1=\left[\begin{array}{cc}16.69& 26.66\\ 27.95& 6.07\\ 12.15& 6.99\\ 7.36& 13.94\end{array}\right]$

非目标情感包括7个样本和2个特征，数据如下：

$w2=\begin{array}{}\left[\begin{array}{cc}-9.88& -24.14\\ 10.58& -6.94\\ 2.22& -13.91\\ -7.93& 3.29\\ -3.88& 5.98\\ -4.86& 1.47\\ 4.28& -1.01\end{array}\right]\end{array}$

Step 2：计算出各类样本均值向量m_i，其中m_i是各类样本的均值，N_i是第ω_i类的样本个数：

$m_i={\left(\frac{1}{N_i}\underset{X\∈W_i}{\mathrm{\Σ}}X\right)}^T$

(16.69+27.95+12.15+7.36)/4＝16.0375；

计算出

$m_1=\left[\begin{array}{c}16.0375\\ 13.4150\end{array}\right]$

$m_2=\left[\begin{array}{c}-1.3529\\ -5.0371\end{array}\right]$

Step 3：计算样本类内离散度矩阵S_i和总类间离散度矩阵S_w：

$S_i=\underset{X\∈{\mathrm{\ω}}_i}{\mathrm{\Σ}}(X-m_i){(X-m_i)}^T,i=\mathrm{1,2}$

S_w＝S₁+S₂

16.69-16.0375＝0.6525；26.66-13.4150＝13.245

$S_1=\left[\begin{array}{cccc}0.6525& 11.9125& -3.8875& -8.6775\\ 13.245& -7.3450& -6.4250& 0.5250\end{array}\right]\×\left[\begin{array}{cc}0.6525& 13.2450\\ 11.9125& -7.345\\ -3.8875& -6.425\\ -8.6775& 0.525\end{array}\right]=\left[\begin{array}{cc}232.7451& -58.4334\\ -58.4334& 270.9353\end{array}\right]$

$S_2=\left[\begin{array}{cc}321.5449& 25.7370\\ 25.7370& 696.5471\end{array}\right]$

$S_w=\left[\begin{array}{cc}554.2900& -32.6965\\ -31.6964& 967.4824\end{array}\right]$

Step 4：计算最佳向量

$w^{*}={\left[\begin{array}{cc}554.2900& -32.6965\\ -32.6965& 967.4824\end{array}\right]}^{-1}\×\left(\left[\begin{array}{c}16.0375+1.3529\\ 13.4150+5.0371\end{array}\right]\right)=\left[\begin{array}{c}0.0326\\ 0.0202\end{array}\right]$

Step 5：计算阈值y₀：

倘若2维样本投影到一维Y空间中，那么它的均值为：

${\tilde{m}}_i=\frac{1}{N_i}\underset{y\∈{\mathrm{\ζ}}_i}{\mathrm{\Σ}}Y=\frac{1}{N_i}\underset{x\∈w_i}{\mathrm{\Σ}}{w}^{T}x=w^T\left(\frac{1}{N_i}\underset{x\∈w_i}{\mathrm{\Σ}}x\right)=w^T{m}_i$

$y0=\frac{{\tilde{m}}_1}{2}+\frac{{\tilde{m}}_2}{2}$

${\tilde{m}}_1=\left[\begin{array}{cc}0.0326& 0.0202\end{array}\right]\×\left[\begin{array}{c}16.0375\\ 13.4150\end{array}\right]=0.7929$

${\tilde{m}}_2=-0.1457$

$y0=\frac{0.7929}{2}+\frac{-0.1457}{2}=0.3236$

Step 6：画出取极大值的向量解如图1所示，：

x＝[-40:0.1:40]为X轴

y＝x*0.0202/0.0326为Y轴

直线为要求出的投影轴，“□”表示目标情感的数据，“×”表示非目标情感的数据。

Step 7：给出测试数据，测试该模型的目标情感和非目标情感的正确识别率ratA和ratB值：

目标情感包括4个样本和2个特征，数据如下：

$p1=\left[\begin{array}{cc}-8.5714& -11.6228\\ 14.1014& 2.5585\\ 7.4654& 8.6988\\ 14.4700& 13.9620\end{array}\right]$

非目标情感包括6个样本和2个特征，数据如下：

$p2=\left[\begin{array}{cc}21.4747& 6.5058\\ 3.9631& 4.7515\\ -2.6728& -9.4298\\ -10.4147& -0.3655\\ -14.2857& -3.8743\\ -14.8387& -6.7982\end{array}\right]$

测试数据分布图如图2所示：

“□”代表阈值，是分割两类情感的，投影大于阈值为目标情感，投影小于阈值为非目标情感，“○”代表目标情感，“☆”代表非目标情感，从图中可看出4个目标情感样本中有1个被识别成了非目标情感，所以

$\mathrm{ratA}=\frac{3}{4}\×100\%=75\%$

6个非目标情感样本中有2个样本被识别成了目标情感，所以

$\mathrm{ratB}=\frac{4}{6}\×100\%=66.67\%$

5.2改进的禁忌搜索算法进行特征选择的步骤如下：

Step 1：将运行“特征选择算法研究综述(毛勇，周晓波.特征选择算法研究综述.模式识别与人工智能[J]，2007，4：V01.20 No.2)”一文中记录的序列后向算法，采用序列后向算法与Fisher分类器结合，进行计算后，将得到的结果中的每个元素用“0”或“1”表示(“0”代表在进行特征选择时该元素没有被选中，“1”代表该元素被选中)，并把该结果放在表L中；j表示表L的行数，设j＝1，为表L的第一行；

Step 2：初始化：禁忌表T＝Ф，设置禁忌长度，设置最大迭代步数；初始解x(j)＝L(j，：)(把表L的第j行作为禁忌搜索算法的第j个空间的初始解)，并把它作为该空间的暂定全局最优解Bestsofar和迭代搜索的起点(当前局部最优解(cand))，计算初始解的适应度函数f(x(j))＝sqrt(ratA²+ratB²)值。

Step 3：单个搜索空间的终止条件：判断是否满足最大迭代步数的停止准则，如满足则终止该空间的计算，并把Bestsofar放在表S的第j行，然后转Step 7，如不满足，则将当前局部最优解(cand)为下次的迭代起点；

Step 4：生成N个候选解集：分别将cand的第i(1≤i≤N，N为特征总维数)位cand(i)值变为1-cand(i)，从第i位的下一位开始查找，遇到值等于1-cand(i)的p位，将其值变为1-cand(p)后结束；若P＝N，则P从1开始继续查找，遇到值等于1-cand(i)的p位，将其值变为1-cand(p)后结束；

Step 5：寻优：计算出每个候选解的适应度函数f的值，从候选解集中选出适应度函数值最大的解，将该解与禁忌表中的解比较，如果该解在禁忌表中，若不满足特赦准则，将较差解作为局部最优解，若满足特赦准则，将禁忌表中的该解提前释放，并作为局部最优解，并将该解的适应度函数值与Bestsofar的适应度函数值比较，如大于，则用该解作为Bestsofar；其中特赦准则指的是，禁忌表中的某个值在迭代过程中作为局部最优解出现规定次数；

如果该解不在禁忌表中，将其中较大的值作为当前局部最优解和Bestsofar的值；

Step 6：更新禁忌表：将当前局部最优解写入禁忌表，转Step 3；

Step 7：转到下一个空间搜索：j＝j+1；如果j不等于特征总维数转Step 2；

Step 8：比较每个空间的解：选出S表中Max(f)对应的解作为最优解(Y)，输出最优解。

下面以原始特征只包含4个特征选择过程为例进行具体介绍：

将运行序列后向算法选出的较好特征组合，放在L表中，j表示L的第j行，此时j＝1：

表4序列后向算法进行特征选择的结果

‘1’表示该特征被选中，‘0’表示该特征没有被选中。

初始化：设置禁忌表为空，空表如表5，设置禁忌长度为2，特赦准则规定的次数为2，设置最大迭代步数为2，将后向算法得出的解“0001”作为暂定的最优解Bestsofar和迭代搜索的起点(当前局部最优解)，采用fisher分类器计算出该初始解的适应度函数值为0.62；

表5禁忌表

空间1的终止条件：判断是否满足最大迭代步数，因为此时未进行迭代，因此不满足；则将当前局部最优解为下次的迭代起点；

根据迭代起点“0001”生成候选解集：

第一个候选集：i＝1，cand(1)的值为0，将cand(1)改为1-0，即1，从第二位开始查找，遇到第4位为1，将其变为0结束；第一个候选集的结果为“1000”；

第二个候选集：i＝2，cand(2)的值为0，将cand(2)改为1-0，即1，从第三位开始查找，遇到第4位为1，将其变为0结束；第二个候选集得结果为“0100”；

第三个候选集：i＝3，cand(3)的值为0，将cand(3)改为1-0，即1，从第四位开始查找，遇到第4位为1，将其变为0结束；第二个候选集得结果为“0010”；

第四个候选集：i＝4，cand(4)的值为1，将cand(4)改为1-1，即0，此时，p的初始值为4，因此重新从第一位开始查找，遇到第1位为0，将其变为1-0，即1结束，；第四个候选集得结果为“1000”；因为第一个候选集与第四个候选集重复，此时只取三个候选集即可。为“1000”，“0100”，“0010”；

寻优：计算候选解集每个解的适应度函数值，“1000”，“0100”，“0010”的适应度函数值分别为0.63，0.55，0.65，可看出“0010”的适应度函数值最大，将“0010”作为当前局部最优解，将“0010”与Bestsofar的适应度函数值比较，发现“0010”适应度函数值大，则将“0010”作为Bestsofar；

更新禁忌表：将“0010”写入禁忌表如下。

表6禁忌表

0	0	1	0

此时只迭代1次，不满足终止条件，将当前局部最优解作为下次迭代的起点，即从“0010”开始进行再次迭代，得到三个候选解集“1000”“0100”“0001”，进行适应度函数值计算后，发现“0010”适应度函数值最大，但是此时“0010”在禁忌表中，选择一个较差的候选解“1000”作为局部最优解，将其与Bestsofar比较，发现“1000”的适应度函数值较小，不替换Bestsofar中的值；

更新禁忌表，将“1000”加入禁忌表，

表7禁忌表

1	0	0	0
				0	0	1	0

对“1000”进行迭代，得到候选解集“0100”“0010”“0001”，进行适应度函数值计算后，发现“0010”适应度函数值最大，将其与禁忌表中的解做比较，发现其已经存在于禁忌表中，满足特赦准则，提前释放“0010”，将其作为局部最优解和迭代的起点，将局部最优解与Bestsofar比较，发现一致，无需替换Bestsofar，并将“0010”再次放入禁忌表，此时，禁忌表更新为

表8禁忌表

0	0	1	0
				1	0	0	0

迭代了3次，符合单个搜索空间的终止条件，把Bestsofar放在表S的第一行，如下表所示：

表9放置各空间最好解的表

转到j＝2的空间搜索：得到第二个空间的解为

表10放置各空间最好解的表

转到j＝3的空间搜索，得到第三个空间的解，

表11放置各空间最好解的表

转到下一个空间搜索：j＝4，j等于特征总维数，结束。

比较每个空间的解：选出“0111”的适应度函数最大，所以将其作为最优解(Y)，输出最优解“0111”。

6.ECG信号的情识别结果

(6-1)相关性分析降维后的特征子集维数

表12相关性降维后的特征子集维数

目标情感	降维后的特征子集维数
		Anger(愤怒)	67
Disgust(厌恶)	65
		Fear(恐惧)	67
Grief(悲伤)	68
		Joy(高兴)	65
Surprise(惊奇)	65

(6-2)ECG信号的情感识别结果

将有效特征或特征组合进行分析比较。随机划分采集得到的数据样本，按照3∶1的比例得到训练、测试。把每类情感生理反应样本分别分为两份：随机抽取其中的3/4构成X_a*n，用于训练分类器；在剩余数据构成T_b*n，用于测试，完成特征选择；其中a、b是样本个数，n是特征维数。并且规定，分类时将情感数据分为目标情感和非目标情感，如恐惧情感的数据为目标情感时，非目标情感即为其它五类(高兴、惊奇、厌恶、悲伤和愤怒)情感的数据。最后总结出能够区分六种情感状态的特征或特征组合，得出六种情感到ECG信号生理特征的映射。

表13相关性降维后基于改进的禁忌搜索算法的ECG信号情感识别

表14基于改进的禁忌搜索算法选择的最优特征子集

击中率(TPR)和虚报率(FPR)来表示，它们分别定义为：

表12是经过相关性分析降维后得到的特征维数，表13经过相关性降维后采用改进的禁忌搜索进行特征选择得到的结果，包括识别六种情感的测试识别率及得到的最优特征子集维数。这些特征能在一定程度上体现情感的变化，对情感的识别有较大的贡献度。表14是经过相关性降维后进行特征选择得到的每种情感的特征名称，实验证明采用表14的特征组合区分某一情感和非某一情感时，识别效果比较好，而且特征维数也较少。因此，找到了识别六种情感状态与心电信号特征的一种映射关系。从以上实验结果可看出ECG包含丰富的情感信息，采用ECG信号用于情感识别研究是可行的，而且经过相关性分析后特征维数减少，运行时间也会减少，得到的实验结果也较好。所以通过相关性分析降维再结合改进的禁忌搜索算法进行特征选择具有明显的优势。

Claims (2)

1.一种用于心电特征选择的禁忌搜索方法，其特征在于：所述方法包括如下步骤：

采用序列后向算法，形成一个N-1行，N列的二维表L，其中N为选择的特征总维数，每一列代表一个特征，每一行称为一个空间，其中第n个空间选择有n个特征，1≤n≤N-1；表内每个元素的值用“0”或“1”表示，“0”代表在进行特征选择时该元素没有被选中，“1”代表该元素被选中；

对每一空间内的选中的特征采用禁忌搜索算法进行求解，得到每个空间的解组成的表S；

选择各个空间中适应度函数最大的作为最终特征选择结果；

其中每个空间中禁忌搜索算法的步骤为：

S1：初始设禁忌表T=Φ，设置禁忌长度，设置最大迭代步数；将此空间内通过序列后向算法得到的值作为初始解,并把它作为该空间的暂定全局最优解Bestsofar和迭代搜索的起点，即当前局部最优解cand，计算初始解的适应度函数

值；

其中ratA和ratB分别代表目标情感的正确识别率和非目标情感的正确识别率；

S2：判断是否满足最大迭代步数的停止准则，如满足则终止该空间的计算，将把Bestsofar放入表S中，结束此空间计算；

如不满足，则将当前局部最优解cand为下次的迭代起点；

S3：生成N个与该空间具有相同特征选择个数的候选解集；

S4：寻优：计算出每个候选解的适应度函数f的值，从候选解集中选出适应度函数值最大的解，

将该解与禁忌表中的解比较；如果该解不在禁忌表中，将其中较大的值作为当前局部最优解和Bestsofar的值；则转S5；

如果该解在禁忌表中，且不满足特赦准则，将适应度函数值第二大的解作为局部最优解；

如果该解在禁忌表中，且满足特赦准则，将禁忌表中的该解提前释放，并作为局部最优解，并将该解的适应度函数值与Bestsofar的适应度函数值比较，如大于，则用该解作为Bestsofar；其中特赦准则指的是，禁忌表中的某个值在迭代过程中作为局部最优解出现规定次数；

S5：更新禁忌表：将当前局部最优解写入禁忌表，转S2。

2.如权利要求1所述的用于心电特征选择的禁忌搜索方法，其特征在于：

所述S3步骤中生成候选解集的方法为：分别将cand的第i位cand(i)值变为1-cand(i)，从第i位的下一位开始查找，遇到值等于1-cand(i)的p位，将其值变为1-cand(p)后结束；若P=N，则P从1开始继续查，遇到值等于1-cand(i)的p位，将其值变为1-cand(p)后结束；

其中1≤i≤N，1≤P≤N,N为特征总维数。

Images