基于MATLAB的脑电信号处理

南京航空航天大学基于Matlab的脑电信号处理姓名陆想想专业领域生物医学工程课程名称数字信号处理二О一三年四月摘要：脑电信号属于非平稳随机信号，且易受到各种噪声干扰。本文基于Matlab仿真系统，主要研究了小波变换在脑电信号处理方面的应用，包括小波变换自动阈值去噪处理、强制去噪处理，以α波为例，提取小波分解得到的各层频率段的信号，并做了一定的分析和评价。关键词：脑电信号；小波变换；去噪重构；频谱分析0引言脑电信号EEG(Electroencephalograph)是人体一种基本生理信号，蕴涵着丰富的生理、心理及病理信息，脑电信号的分析及处理无论是在临床上对一些脑疾病的诊断和治疗，还是在脑认知科学研究领域都是十分重要的。由于脑电信号的非平稳性且极易受到各种噪声干扰，特别是工频干扰。因此消除原始脑电数据中的噪声，更好地获取反映大脑活动和状态的有用信息是进行脑电分析的一个重要前提。本文的研究目的是利用脑电采集仪器获得的脑电信号，利用Fourier变换、小波变换等方法对脑电信号进行分析处理，以提取脑电信号α波的“梭形”节律，并对脑电信号进行功率谱分析和去噪重构。1实验原理和方法1.1实验原理1.1.1脑电信号根据频率和振幅的不同，可以将脑电波分为4种基本类型[1]，即δ波、θ波、α波、β波。4种波形的起源和功能也不相同，如图1所示。图1脑电图的四种基本波形α波的频率为8~13Hz，振幅为为20~100µV，它是节律性脑电波中最明显的波，整个皮层均可产生α波。正常成人在清醒、安静、闭目时，波幅呈现有小变大，再由大变小，如此反复进行，形成所谓α节律的“梭形”。每一“梭形”持续时间约为1~2s。当被试者睁眼、警觉、思考问题或接受其他刺激时，α波立即消失而代之以快波，这种现象称之为“α波阻断”。一般认为，α波是大脑皮质处于清醒安静状态时电活动的主要表现。β波的频率是18~30Hz，振幅为5~20µV，是一种快波。β波的出现以为着大脑比较兴奋。θ波的频率是4~7Hz，振幅为10~50µV，它是在困倦时，神经系统处于抑制状态时所记录的波形。δ波在睡眠、深度麻醉、缺氧或大脑有器质性病变时出现，频率是1~3.5Hz，振幅为20~200µV。1.1.2小波变换小波变换的概念是由从事石油信号处理的法国工程师J.Morlet在1974年首先提出的。与Fourier变换相比，小波变换是时间（空间）频率的局部化分析，它通过伸缩平移运算，对信号逐步进行多尺度细化，最终达到高频处时间细分，低频处频率细分，能自动适应时频信号分析的要求，从而可以聚焦到信号的任意细节，解决了Fourier变换的困难问题。在噪声中如何准确地检测到信号一直是信号处理领域所关心的内容。小波变换是一种信号的时间一尺度分析方法，由于具有多分辨率分析的特点，良好的时频局部化特性，能够对各种时变信号进行有效的分解，从而较好地将信号与噪声加以分离，获得满意的去噪效果。小波分析去噪处理的方法一般有三种：默认阈值去噪处理。该方法是利用ddencmp函数生成信号的默认阈值，然后利用wdencmp函数进行去噪处理。给定阈值进行去噪处理。在实际的去噪过程中，阈值往往可以通过经验公式获得，并且这种阈值比默认阈值可信度高。在进行量化处理时可利用wthresh函数。强制去噪处理。该方法是将小波分解结构中的高频或者低频系数设置为0，即滤掉所有高频部分或低频部分。这种方法比较简单，且去噪得到的信号比较平滑，但是容易丢失信号中的有用成分[2]。本文采用了两种去噪方法，并分析比较了他们的去噪效果。1.2实验方法与步骤1.2.1脑电信号的读取本文使用的脑电数据是使用南航生物医学光子学实验室的脑电采集系统采集获得的，原始数据格式为.eeg。为了方便在Matlab环境下对数据进行分析，将文件转换为.txt格式。脑电采集使用的是16通道，采样频率为256Hz，文件中存储的数据的形式为数据点数×通道数。实验中选取了第14通道的前8000个数据点作为样本进行分析。由于采样时间是256Hz，所以这段信号的持续时间大约是32秒。1.2.2信号的频域和功率谱分析为了研究脑电信号中不同频率信号的能量分布以及变化情况，首先对样本信号进行Fourier变换，得到频域图。然后进一步对信号做功率谱分析，得到功率谱图，从功率谱图中，可以直观的观察到不同频率信号的能量分布情况。由于脑电数据是在被采集者安静清醒的状态下采集得到的，理论上α波应该占主导地位。1.2.3信号的小波变换及重构基于小波变换降噪处理的方法通常有3个步骤：首先是将信号进行n层小波变换，得到小波系数；然后在小波变换域上利用信号与噪声的不同特性，对小波变换进行阈值化处理，把噪声从信号中区分开来(主要是对高频系数进行阈值化处理；最后是利用重构算法重构信号。小波变换去噪的效果主要取决于对含噪信号的噪声估计方法以及所采用的小波函数[4]。本文使用其中两种去噪方法。第一种是默认阈值去噪，先对样本信号进行8层小波分解，使用的小波函数是db8，然后利用ddencmp函数生成信号的默认阈值，最后利用wdencmp函数进行去噪处理，得到去噪的信号。另一种是强制去噪，同样地，先对样本信号进行8层小波分解，采用的小波函数也是db8。然后提取各层小波系数，再将脑电信号频率以外的几层小波系数置零，得到重构的各层小波系数。最后由重构的小波系数得到重构的信号。该信号中除去了某些频率的信号，起到了去噪的效果。2仿真结果及分析2.1原始脑电数据的读取和显示采集到的脑电信号文件为data.txt，调用eeg_load.m文件，即可绘制出脑电样本信号图，如下图2所示。图2脑电样本信号图2.2脑电信号频谱图及功率谱图的绘制首先调用eeg_fft.m文件，原理是对样本信号进行傅立叶变换[3]，即可获得样本信号的频谱图，如下图3所示。图3信号频谱图从频谱图中可以看出，低频信号和11Hz左右的信号特别强，25Hz以上的信号几乎为零。由于α波的频率为8~13Hz，由此可以知道，该信号中α波比较多。在脑电采集过程中，被采集者没有处于深度睡眠状态，因此接近0Hz的低频信号可以确定为噪声。15~25Hz频段的信号很微弱，因此可以以判断出，信号中几乎没有β波。谱估计是数字信号处理的主要内容之一，主要研究信号在频域中的各种特征，目的是根据有限数据在频域内提取被淹没在噪声中的有用信号。调用eeg_psd,m文件，可以绘制信号的功率谱图，如下图4所示。图4信号功率谱图原始的功率谱图是关于f=128Hz对称的，为了便于分析，截取了0~65Hz的一段。信号的功率谱显示，在10Hz作用频率处能量出现峰值，表明在平静状态下采集的脑电信号中α波能量最大，符合生理学的研究结论。2.3脑电信号小波分解各层重构波形实验选取db8小波对前述采集的样本信号进行分解，调用wavelet_dec.m文件，分解得到的各层信号如下图5所示。图5脑电信号经小波分解后的各层分量由于采样频率是256Hz，所以对于自带的频率范围是：]2,2[],2,2[],2,2[],2,2[],2,2[],2,2[],2,2[],2,0[22334455667788sssssssssssssssfffffffffffffff相应的各子带的频率成分如下表所示：子带信号频率范围主要信号成分a80Hz~1Hz低频干扰d81Hz~2Hzδ波d72Hz~4Hδ波d64Hz~8Hzθ波d58Hz~16Hzα波d416Hz~32Hzβ波d332Hz~64Hz高频噪声d264Hz~128Hz高频噪声d1128Hz~256Hz高频噪声表1小波分解后脑电信号子带的频率范围一般情况下，a8子带内是低频干扰，d3、d2、d1子带内是高频噪声，d8、d7、d6、d5、d4子带内是脑电信号，但也可能混有一定的噪声，需要根据实际情况来分析。2.4脑电信号节律提取α波是节律性脑电波中最明显的波。由以上小波分解得到的各层分量分析可知，α波主要集中在d5子带内。所以选取了d5子带内第5500~8000之间的2500个数据点，将信号图绘制出来，如下图6所示。图6提取出来的α波节律从上图中可以明显看出α波波幅由小到大，再由大到小作规律性变化，呈棱状图形。2.5小波分解去噪和重构波形2.5.1小波变换默认阈值去噪处理调用noise_reduce.m文件，可以实现脑电信号的小波变换默认阈值去噪处理，原始数据及去噪处理结果对比如下图7所示。图7原始信号与小波默认阈值去噪结果图的对比从原始脑电信号与去噪处理后的效果来看，经去噪处理后的信号高频信号有所减少。2.5.2小波变换强制去噪处理调用wavelet_rec.m文件，绘制小波变换强制去噪处理之后的信号如下图8所示。图8原始信号与小波强制去噪结果图的对比该图表示的是除去了低频干扰和高频噪声之后的结果，从与原始信号的对比中可以看出，高频噪声很明显被消除了。但由于考虑到实际情况，16~32Hz子带内极少有β波，大部分为噪声，所以把这个频带内的信号也全部清除了。3讨论基于FFT变换的功率谱估计适用于平稳时间信号分析,计算结果只能反映信号段总体的平均功率分布情况,不包含信号的任何时域变化信息,并且谱估计的频率分辨率与所采用的信号长度成正比，即受数据点数目的影响。基于小波变换的去噪方法，对非平稳信号去噪，要比传统的滤波去噪声得到的效果好，主要是由于传统的滤波器都具有低通性，对需要分析在每个时刻含有不同频率成分的非平稳信号来说，是很难对它进行匹配分析。而小波变换具有多分辨率，且在时频域都具有局部性，因此很适合用来分析非平稳信号。在用小波分析来进行去噪的关键在于阈值的选取，如果阈值选取的太高，会使信号失去太多细节，使信号失真，如果阈值选得太低，又不能达到去噪的目的。4结论本文利用实测的原始脑电信号,对脑电信号的处理方法与结果进行了一定的分析和评价,以期为脑电信号处理及特征提取提供一定的理论参考和分析依据。目前人们也尝试用非线性处理方法、神经网络的方法、时频结合等等现代的方法来处理脑电信号,相信这些方法会为脑认知以及医学的发展作出贡献[5]。脑电信号属于非平稳随机信号，小波分析的方法可以直接对信号的某些频率分量进行观察或者提取出有用的特征信号，为脑电信号的测量与分析提供了非常好的前景。参考文献:[1]余学飞.现代医学电子仪器原理与设计[M].广州:华南理工大学出版社,2007.133~134.[2]张德峰.MATLAB小波分析[M].北京:机械工业出版,2012.[3]万永革.数字信号处理的MATLAB实现[M].北京:科学出版社,2012.[4]于兰兰.基于小波变换的脑电信号去噪处理[J].南昌大学学报(理科版),2007,31:75~77.[5]谢松云,张振中,杨金孝,张坤.脑电信号的若干处理方法研究与评价[J].计算机仿真,2007,24(2):326~330.