毕业设计-脑电信号预处理方法研究-马孝龙

本科毕业设计（论文）

脑电信号预处理方法研究

马孝龙

燕山大学

2013 年 6 月

本科毕业设计（论文）

脑电信号预处理方法研究

学院：信息科学与工程学院

专业：通信工程

学生姓名：马孝龙

学号： 0901******** 指导教师：李段

答辩日期： 2013-6-24

燕山大学毕业设计（论文）任务书

摘要

摘要

脑电信号(EEG)是脑神经细胞电生理活动在大脑皮层或头皮表面的总体反映，其中包含了大量的生理和病理信息,并可以用许多特征量来描述其特征信号。通过脑电分析来认识脑的活动是一种有效的无创手段。但是，脑电信号中常常混杂有多种噪声，为了便于阅读和分析脑电信号，必须对脑电信号进行去噪处理。

本文针对脑电信号中的眼电伪迹去除的问题，采用了小波阈值法对脑电信号进行去噪处理。以小波阈值降噪为基础，本文采用db5小波对脑电信进行5尺度分解，其次在对应的子带上采用软阈值法，且阈值原则为对数长度阈值法对小波进行阈值去噪，然后对去噪后的信号进行重构，最后比较原始信号和去噪后的信号。实验结果表明，用小波阈值法可以有效去除脑电信号中的眼电干扰。

关键词脑电信号; 小波变换; 去噪

Abstract

Abstract

The Electroencephalograph (EEG) is the electro physiological activity of brain cells, it is reflected in the overall pallium or scalp surface, which contains a large number of physiological and pathological information, and can describe the characteristic signal with many features. The EEG analysis to understanding brain activity is an effective noninvasive method. However, EEG is often mixed with a variety of noise, in order to facilitate the reading and analysis of EEG signal denoising, must carry on eeg.

According to the problem that removing the EEG of ocular artifacts, this article uses the wavelet threshold denoising method to remove the eye electrical noise (EOG) interference in the EEG. Based on wavelet threshold denoising, firstly, this article selects the db5 wavelet scale decomposition of EEG of 5. Secondly, it uses soft threshold method in the corresponding subband. Then the soft threshold and wavelet reconstruction algorithm for denoising, and the threshold principle for the fixed threshold method of wavelet threshold denoising. Then reconstruct the signal after denoising. Finally, comparison of the original signal and the denoised signal. From the result of experiments, we can see that it is perfectly to remove the eye electrical noise (EOG) by using the wavelet threshold.

Keywords EEG; Wavelet transform; noise rejection

目录

摘要 ....................................................................................................................... I Abstract ................................................................................................................ II 第1章绪论 .. (1)

1.1课题背景 (1)

1.2国内外研究现状与成果 (2)

1.3本文研究内容与结构安排 (4)

第2章小波变换在信号处理中的应用 (5)

2.1小波变换 (5)

2.1.1连续小波变换 (5)

2.1.2离散小波变换 (7)

2.1.3二进小波变换 (8)

2.1.4常用小波 (9)

2.2小波去噪 (12)

2.2.1尺度间相关性的去噪方法 (13)

2.2.2阈值去噪法 (14)

2.3MATLAB小波在信号去噪中的应用 (16)

2.3.1信号分解 (17)

2.3.2信号重构 (18)

2.3.3信号阈值去噪 (20)

2.4本章小结 (21)

第3章眼电噪声去除方法研究 (23)

3.1脑电信号与眼电干扰简介 (23)

3.2小波阈值去噪原理 (25)

3.3参数选择 (26)

3.3.1小波基的选择 (26)

3.3.2 小波分解层数的选择 (27)

3.4实验结果仿真与分析 (27)

3.5本章小结 (32)

结论 (33)

参考文献 (35)

致谢 (37)

附录1 (39)

附录2 (44)

附录3 (49)

附录4 (54)

第1章绪论

第1章绪论

1.1 课题背景

大脑是一个非线性、非平稳的，且由大量神经元群相互连接形成的复杂动力学系统，人们获取大脑信息的手段基本上分为三种：通过脑结构的解剖图像获得静态的信息，例如广泛使用的电子计算机X射线断层扫描技术和磁共振成像；通过功能活动时的代谢活动图像而获得动态的信息，例如功能磁共振成像技术和正电子发射断层扫描；通过大脑对信息处理的反应而获得脑的动态信息，这类手段有局部场电位，脑电图和脑磁图。

其中，EEG 是大脑组织中大量的神经元群突触后电流的在大脑皮层或头皮表面的综合表现，EEG 具有简易、无创、高时间分辨率等优点，脑电信号处理可以获取大量的功能与疾病信息，例如通过EEG检查，可以对癫痫、脑肿瘤、颅脑损伤、中枢神经系统疾病，以及精神疾病和智能障碍等疾病提供信息，从而对脑功能分析及疾病诊断、预后和治疗提供有效的方法；也可利用人对不同的感觉、运动或认知活动的脑电的不同，通过对脑电信号的有效的提取和分类达到某种控制目的。

10个神经元组成，典型的神经大脑皮层是大脑信息处理的中枢，约由11

元由细胞体、树突和轴突三部分组成，当神经元接受的其它神经元发来的信号的和超过阈值时，神经元就会发出神经脉冲，称作动作电位，神经元的轴突终末与其它神经元的接受表面形成各种突触，实现神经元间的信息交换，一个典型的神经元可能有1000-10000个突触，能接受来自1000个其它神经元的信息，突触后电位在其产生部位的附近按空间总和规律叠加起来，脑电信号就是由皮层内大量神经元突触后电位(包括兴奋性突触后电位和抑制性突触后电位两种)总和所形成的，是许多神经元共同活动的结果。

随着信息科学的不断深入，脑科学的研究取得了新的进展。能够很好反应脑信息的脑电(electroencephalogram, EEG)是大脑神经元突触后电位的综合结果，是大脑活动产生电场经容积导体传导后在脑皮上的电位分布。现代的脑科学研究表明，多导的脑电信号蕴含着多种生理现象且能反映更多的动态信息。多导脑电信号就是从头皮上按照一定标准放置的多个电极处同步采

燕山大学本科生毕业设计（论文）

集的电信号。因为脑电是通过安放在脑皮层不同位置的电极测量得到的，而电极之间的距离只有几厘米，所以这些电极不仅记录了脑内神经元电活动，同时也记录了各种干扰信号，如工频干扰、眨眼、眼球运动、心电干扰、呼吸波干扰和肌电干扰等等[1,2]。因此，为了得到真实的脑电数据，使之能够在临床上合理正确地做出判断，就必须对脑电测量数据进行预处理。

1.2 国内外研究现状与成果

眼电信号是脑电信号的主要干扰源，它们能够改变脑电信号的电场分靠使获得的脑电信号不能真实的反映大脑的活动。眼电伪迹去除的基本原则是在保证脑电信号完整的前提下尽可能的去除眼电伪迹。这个基本原则是评价各种眼电伪迹去除方法的重要依据。目前，常用脑电信号去噪方法有独立分量分析(ICA)的脑电信号去噪方法,有主成分分析算法的脑电信号去噪方法，基于自适应滤波的脑电信号去噪方法,基于陷波器的脑电信号去噪方法,基于小波变换的脑电信号降噪方法等多种方法。

(1) 基于独立分量分析(ICA)的脑电信号去噪方法[3,4]工频噪声、心电伪迹以及脑电波源信号之间的关系是统计独立的,满足ICA方法的分离条件,可将脑电信号去噪问题转化为独立分量分离问题,通过构造与工频噪声频率相同的正交正弦和余弦信号作为对工频噪声的参考信号,将构造的两个参考信号和心电信号以及含噪脑电信号作为ICA中混合矩阵的输入信号,采用收敛速度快的Fast ICA算法把脑电信号中的工频噪声和心电伪迹作为独立信号分离出去,得到去噪后的脑电信号。结果通过ICA方法对噪声进行分离后,脑电信号中的两种噪声基本被消除,并且可很好的保留脑电信号有用成分。

(2) 基于主成分分析的脑电信号去噪方法此算法去除眼电伪迹的基本思想是：假设脑电信号与眼电信号彼此正交，将原始数据投影到彼此正交的特征空间，从而实现原始信号的特征分离，识别并去除眼电伪迹，逆向投影恢复原始数据。

(3) 基于陷波器的脑电信号去噪方法[5,6]杨帮华等人共同提出基于陷波器和小波变换去除自发脑电信号噪声的方法。该方法是针对自发脑电信号的特点构造了一种零通小波，根据变换尺度参数与频率的对应关系，有选择

第1章绪论

地某些感兴趣的尺度信号以去除噪声。这种方法可以很好地去除工频干扰外的噪声信号，因此，进行小波去噪后，再结合陷波滤波器去除50Hz的工频干扰，这两种方法结合可以有效地去除干扰。

(4) 基于自适应滤波的脑电信号去噪方法早在20世纪40年代，就对平稳随机信号建立了维纳滤波理论。根据有用信号和干扰噪声的统计特性，以线性最小均方误差估计准则所涉及的最佳滤波器，成为维纳滤波器。这种滤波器能最大限度的滤除干扰噪声，提取有用信号，然而，对随机非平稳的信号，它的实际应用受到了限制。1967年WdrowB.等提出了自适应滤波理论，可使自适应滤波系统的参数自动的调整而达到最佳状况，而且在设计时，只需要很少的或根本不需要任何关于信号与噪声的先验统计知识，这种滤波器的实现差不多像维纳滤波器那样简单，而滤波性能几乎如卡尔曼滤波器一样好。

考虑到脑电信号源的随机非平稳性，用自适应滤波方法去除脑电噪声很合适。到目前为止已有很多自适应滤波对脑电数据处理的研究成果，其中李翠芳等人共同提出基于小波分析的EEG信号自适应去噪的应用[7]。根据短时动态信号与平稳背景噪声的特征区别，对输入信号进行白化预处理，以时间序列的AR模型理论为依据，导出背景噪声白化滤波器的结构。将小波变换与自适应滤波器相结合，对经白化处理后的信号进行自适应去噪。

(5) 基于小波变换的脑电信号降噪方法[8,9] 近年来，随着电子技术的迅猛发展，信息获取的手段、进度、速度都有很大的提高。特别是非平稳信号分析理论上的一系列重大进展为非平稳信号提出了新的处理与分析手段。小波分析理论则是这一系列重大进展中的一个。它不仅继承和发展了窗口傅里叶变换的局部化思想,而且克服了窗口大小不随频率变化, 缺乏离散正交基的缺点。小波变换在时间和频率上都有很好的局部性，在分析低频信号时，其时间窗很大，而在分析高频信号时，时间窗很小，这使得小波变换非常适合于时频分析。小波变换对于信号的高频成分使用逐渐尖锐的时间分辨率以便移近观察信号的快变成分，对于低频成分使用逐渐尖锐的频率分辨率以便移近观察信号的慢变成分（整体变化趋势）。借助时频局部分析性，小波分析理论已经成为信息去噪处理的一个重要工具。

燕山大学本科生毕业设计（论文）

目前，小波分析消噪的方法[2,5]大概可以分为三大类：第一类方法是对含噪信号做小波变换之后，计算相邻尺度间小波系数的相关性，根据相关性的大小区别小波系数的类型，从而进行取舍，然后直接重构信号；第二类方法是基于小波变换模极大值原理，最初由Mallat提出，即根据信号和噪声在小波变换各尺度上的不同传播特性，剔除由噪声产生的模极大值点，保留信号所对应的模极大值点，然后利用余下的模极大值点重构小波系数，进而恢复信号；第三类方法是Donoh提出的阈值方法。该方法认为信号对应的小波系数包含有信号的重要信息，幅值较大，但数目较少；噪声对应的小波系数则是均匀分布的，个数较多，但幅值较小。基于此思想，把绝对值较小的系数置零，而让绝对值较大的系数保留或者收缩(分别对应于硬闽值法和软阈值法)，得到估计小波系数。此方法在Besov空间上可以得到最佳估计值，所以成为目前应用最为广泛的小波去噪方法之一。

小波分析属于时频分析的一种，是从傅里叶分析中发展而来的，但是又优于傅里叶分析。而作为傅里叶分析思想发展的小波分析方法是一种窗口大小，即窗口面积固定、但窗口的形状可变、时间窗和频率窗都可改变的时频局部化分析方法，即在低频部分具有较高的频率分辨率和较低的时间分辨率，在高频部分具有较高的时间分辨率和较低的频率分辨率，很适于探测正常信号中夹带的瞬态反常现象并展示其成分，在时频域都有表征信号局部特征的能力，具有多分辨分析的特点。

1.3本文研究内容与结构安排

本文首先介绍了脑电信号的背景知识，然后介绍了小波理论并利用阈值法滤除脑电信号中的眼电干扰，通过编程成功完成信号分解、去噪、重构等内容。

论文第1章为绪论，介绍了脑电信号的特点、性质、分类，小波变换的由来，以及本课题研究的主要内容和意义；第2章就小波变换在信号处理方面的知识加以介绍，为课题后续研究打下基础。第3章分析了已有方法，比较优劣，并提出了基于小波变换的脑电信号预处理研究方法，以及对脑电信号中夹杂的眼电信号的阈值法去噪。

第2章 小波变换在信号处理中的应用

第2章 小波变换在信号处理中的应用

小波变换的概念是1984年法国地球物理学家J. Morlet 在分析处理地球物理勘探资料时提出的。小波变换的数学基础的19世纪的傅立叶变换，其后理论物理学家 A. Grossman 采用平移和伸缩不变性建立了小波变换的理论体系。1988年比利时数学家I. Mallat 提出了多分辨率分析的概念，统一了在此之前的各种构造小波的方法，特别是提出了二进小波变换的快速算法，使得小波变换完全走向实用化。

小波变换采用改变时间-频率窗口形状的方法，很好的解决了时间分辨率和频率分辨率的矛盾，在时间域和频率域都具有很好的局部化性质。对信号中的低频成分采用宽的时间窗，得到高的频率分辨率；对信号中的高频成分采用窄的时间窗，得到低的频率分辨率。小波变换的这种适应性，使它在工程技术和信号处理方面得到了广泛的应用。

2.1 小波变换

2.1.1 连续小波变换

设函数)(∈)(2R L t ψ,满足下述条件

?=R dt t 0)(ψ (2-1)

称)(Ψt 为基本小波。引入尺度因子a 和平移因子b, a 和b 满足：

0,≠∈a R b a 且

将基本小波进行伸缩和平移，得到下列函数族

)(||)(2/1,a

b t a t b a -ψ==ψ (2-2) 称)(,t b a ψ为分析小波。系数a /1为归一化常数，它使得所有对尺度a 和位移参数b ，下式成立

t d t dt t R R b a b a ??ψψ=ψ22,2,|)(||)(||||| (2-3)

通常取

?=ψR dt t 1|)(|2，其意义的)(t ψ具有单位向量。

函数)()(2R L t f =的连续小波变换(CWT)定义为 1/2,,,()()()()a b a b a b R R t b CWT f t t dt a

f t t a -=ψ=ψ?? (2-4)

燕山大学本科生毕业设计（论文）

若基本小波)(Ψt 满足下述条件：

∞)(2

<ψ?R d ωωω (2-5)

则连续小波变换CWT a,b 存在逆变换，逆变换公式为

2b a,,01-a dadb CWT )(1)(?ψ=?∞ψt a

C t f b a (2-6) ?ψ=ψR d ωωω2

)(C

(2-7) 称式（2-5）为容许条件，称满足容许条件的小波为容许小波。

(1)，线性 连续小波变换是线性变换，即一个函数的连续小波变换等于该函数的分量的变换和。公式如下：

)(+)(=)(21t f t f t f

CWT t f ?)(

11)(CWT t f ?

22?)(CWT t f

则 CWT=CWT 1+WT 2

(2)，时移性

b a CWT t f ,?)(

则 0t -,0?)t -(b a C W T t f

(3)，时标定理

b a CWT t f ,?)(

cb ca CWT t f ,?)c (

(4)，微分运算

dt t t

t f t t f b a m m m m m )( ??)()1-() ?)(?(CWT ,b ,a ψ=?∞∞- (2-8)

(5)，能量守恒

???∞∞∞∞∞ψ=-0-22a 1-2a dadb |CWT |)(b C t f ， (2-9)

(6)，冗余度 连续小波变换是把一维信号变换到二维空间b a CWT t f ,?)(，因此小波变换中存在多余信息，成为冗余度(Redundancy)。因而小波变换公式不是唯一

第2章 小波变换在信号处理中的应用

的。从分析小波)(

||)(2/1,a

a t

b a ψ==ψ角度看，)(,t b a ψ是一族超完备基函数，他们之间是线性相关的。度向量冗余度成为再生核1212(,,,)K a a b b ， dt a b a a b a ),,,(2

21/2-2111/2-1-1-2121-ψ-ψ=?∞∞ψt t C b b a a K (2-10) 再生核就是小波本省的变换。

再生核),,,(2121b b a a K 度量了小波变换二维空间里两点）（11,a b 与）（22,a b 之间的相关性大小。再生核K 作用于小波变换CWT a,b ,仍得到CWT a,b

??∞∞-∞=02b',a'112/1-b a,'

''CWT )',,',(CWT 'a da db b b a a K a (2-11) 2.1.2 离散小波变换

信号)(t f 的连续小波变换为

dt t t f CWT b a R b a )()(,,∫ψ=

对尺度因子a 和平移参数b 进行如下的离散采样：

m a a 0=

Z m a ∈>,00 m a nb b 00=

Z n b ∈,0> 则小波)(Ψ,t b a 变为 )n ()(0-m 02/-m 0,b t a a t n m -=ψ

离散小波变换定义为 ?ψ=R n m b a dt t t f DWT )()(,, (2-12)

写成内积的形式有

)(),(,,t t f DWT n m b a ψ= Z n m ∈,

(2-13)

对于离散小波做如下三点说明： 等Q 结构离散化 由于基本小波的等Q 性质，对参数也做等Q 结构离散

化，即a 增大时，a 的间隔也增大，所以取m a a 0=，a 延时时间也增大，故

取b 为m a nb 00。参数离散化的小波为

)n ()(0-m 0-m/20

,b t a a t n m -=ψ (2-14) 离散小波变换实际上仍然是一系列带通滤波器，只是带通滤波器的中心频率也带宽由于a 的离散采样率而为一系列的离散值。但仍然保持恒Q 性质。滤波后的输出也因ｂ的离散采样率而为离散值。

燕山大学本科生毕业设计（论文）

离散小波的重构。 根据前述框架理论结果，当{n m ,ψ}为)(2R L 框架时，可由离散小波变换

n m DWT ,恢复出原信号f(t),其重构公式为 )()(,~,,∑t DWT t f n m n m n m ψ=

(2-15)

)(,~

t n m ψ为)(,t n m ψ的对偶框架，而 dt t t f DWT n m R b a )()(,,ψ=?

(2-16) 离散小波的冗余性。离散小波框架{)(,t n m ψ}存在冗余性，所以离散小波变换的冗余变换。但是，当框架界A=B=1时，{)(,t n m ψ}就成了)(2R L 的正交基。这时，信号f(t)的离散小波变换就相当于正交分解，称为正交离散小波变换，它是无冗余变换。

2.1.3 二进小波变换

在离散小波变换中，一种方便的离散方法是取a 0=2，所得到的小波和小

波变换称为二进小波和二进小波变换。如果再去b 0=1.称其为二进正交小波

和二进正交小波变换。

设)(∈)(2R L t ψ，若存在常数A 和B ，∞<≤

B A k k

≤ψ<∑2

)2(ω (2-17) 则称)(t ψ为二进小波。条件(2-17)成为稳定条件。若A=B ，则称为最稳定条件。二进小波是容许小波。

对于二进小波{)2(2)(k --k 2k t t ψ=ψ}、)()(2R L t f ∈的二进正交小波

变换定义为 dt t t f k R

)2b ()(2DWT k -2k -ψ=? (2-18) 很容易得到k 2DWT 的傅立叶变换为

)2()(DWT k 2ωωk F ψ= (2-19)

利用式(2-17)和(2-18)有：

2

2k 22)F(B ||)(||)F(A ωωω≤≤∑k DWT (2-20)

式(2-19)表明二进小波的框架。根据框架理论结果，二进小波变换k 2DWT 可

第2章 小波变换在信号处理中的应用

重构为f(t)。由框架重构公式知道，需要再给出重构小波，为此定义下述方程：

1)2()2

(=ψ∑k k k ωχω (2-21)

式中)(t χ为重构小波，通常取)(ωχ为

2|)2

(|)

()(∑ψψ=k k ωωωχ

因此下面给出二进小波的重构公式

?∑?

?????-?=R k db t b DWT t f k )]b (2[)()(k -2χ (2-22) 二进正交小波变换

设)(∈)(2R L t ψ，且满足

1|)2(|2=ψ∑k k ω

(2-23)

)(t ψ为正交二进小波。尺度因子和平移参数按二进制离散，1,200==b a ，二进正交小波变换为

)n 2(2)(-2/-k ,-ψ=ψt t k n m (2-24)

2.1.4 常用小波

小波分析的基（小波函数）不是唯一存在的，所有满足小波条件的函数都可以作为小波函数，那么小波函数的选取就成了十分重要的问题，实际选取小波的标准主要有一下三种：

(1) 自相似原则：如果选择的小波对信号有一定相似性，则变换后的能量就比较集中，可以有效减少计算量。

(2) 判别函数：针对某类问题，找出一些关键性的技术指标，得到一个判别函数，将各种小波函数代入其中，得到一个最优准则。

(3) 支集长度：大部分应选择支集长度在5~9之间的小波，因为支集长度太长会产生边界问题，支集太短则不利于能量集中。

下面介绍几种常用的小波：

Harr 小波是小波分析发展过程中用的最早的小波，也是最简单的小波（如图2-1所示），Haar 小波是一组正交归一的函数集，是一个时域不连续

燕山大学本科生毕业设计（论文）

的二进制小波函数，其解析式表示如下：

)252(,012/1,

12/10,1)(-?????≤≤-≤≤=ψelse x x x Daubechies 小波是由著名的小波学者Ingrid Daubechies 创造的，她发明的紧支集小波是小波领域的里程碑，使得小波的研究由理论转入实践。

Daubechies 小波简写为dbN ，其中N 表示阶数，db 是小波名字的前缀，它在时域上具有有限支撑。支集长度和滤波器长度都是2N 左右，消失矩为N ，可见这个序列的小波扩展性比较好，可以比较灵活的权衡增加支集长度带来的边界问题。函数图形如图2-2。

Coiflet 小波 Daubechies 提出的另一个小波系，有更长的直接支集长度（6N-1）和更大的消失矩阵(2N)，对称性比较好。其函数图形如图2-3.

Morlet 小波 Morlet 小波是一个具有解析表达式的小波 (如图2-4)，但它不具备正交性，所以只能满足连续小波的可允许条件，不存在紧支集，不能做离散小波变换和正交小波变换。其解析表达式为：

)5c o s ()(/2-x 2

x Ce x =ψ （2-26）

墨西哥帽（Mexican Hat ）小波 墨西哥帽小波又称为Marr 小波(如图2-5)，类似于Morlet 小波。它同样具有解析表达式，但不存在尺度函数，不具有正交性。它是高斯函数的二阶导数（加负号），在时域和频域都有很好的局部化特性。其解析表达式为：

/2x -24/1-2)x 1(32)(e x -=πψ （2-27）

图2-1 Harr 小波函数图

第2章小波变换在信号处理中的应用

图2-2 db5小波函数

图2-3 coif3小波函数

图2-4 morl小波函数

燕山大学本科生毕业设计（论文）

图2-5 墨西哥草帽小波

symlets morlet 小波函数Harr Daubechie

s

正交性有有有无

双正交性有有有无

紧支撑性有有有无

支撑长度 1 2N-1 2N-1 有限长度对称性对称近似对称近似对称对称

消失矩阶数 1 N N 无

2.2 小波去噪

一般来说，实际采集到的信号都带有噪声，在对信号进行进一步处理之前，需要将有效的信号提取出来。由于噪声的存在，有用信号常常被淹没在强大的噪声之中，因此，去噪和滤波就成了信号处理中的一个重要的内容。目前，人们根据噪声的统计特性和频谱分布规律，开发了多种多样的信号去噪方法。小波去噪法就是其中的一大类，且小波变换具有低熵性、多分辨率性质、去相关性、小波基选择多样性，这为其取得广泛应用打下了坚实的基础。信号去噪的准则有两个：

(1) 光滑性在大部分情况下，降噪后的信号应该至少和原信号具有同等的光滑性。

(2) 相似性降噪后的信号和原信号的方差估计应该是最坏情况下的最小值。