基于脑电的平移不变量小波去噪

1、时频分析与小波变换专题训练报告题 目平移不变量小波去噪方法学 院自动化学院专 业控制理论与控制工程学 号131060024学生姓名马秉宇指导教师佘青山1.前言小波变换的低熵性、多分辨率性、去相关性和小波选择的多样系，使其广泛应用于含噪信号消噪领域。1994年，Donoho等提出了包括软阈值法和硬阈值法的小波阈值收缩法(WaveShrink)方法，该方法在最小均方误差意义下可达近似最优，并且可取得较好的消噪效果。软阈值法消噪得到的信号整体连续性好，不会产生附和的振荡，但是与原信号的逼近程度差；硬阈值法消噪得到的信号与原信号逼近程度好，但会产生附加的振荡，不论采用软阈值方法还是硬阈值方法消噪，在有

2、些情况下，阈值法去噪后信号的某些不连续点附近和信号的快速变化点处，会出现伪吉布斯现象，即在这些点去噪信号会在一个特定的目标水平上下跳变，在一定程度影响了消噪效果。平移不变量(translation-invariant,TI)小波消噪可抑制伪吉布斯现象，如对含噪的Blocks、Bumps、HeaviSine和Doppler信号进行平移不变量小波消噪处理，较好地抑制了伪吉布斯现象；李肃义等对心电信号进行平移不变量消噪处理，达到了抑制伪吉布斯现象的效果，为了提高消噪后信号的连续性和原信号的逼近度，综合软、硬阈值小波消噪信号处理的优势，同时使处理后的信号不再出现伪吉布斯现象，本文将这种方法应用于脑电信

3、号的消噪处理，实验结果表明，该方法可以有效地提高信噪比，降低均方误差，并能较好地保留EEG信号有效特征。2.平移不变量小波去噪原理2.1 正交小波变换1998年，Mallet在构造正交小波基时提出了多分辨率分析的概念，并给出了正交小波变换快速算法，即Mallet算法，根据过分辨率分析的理论，若为信号的离散采样数据，则信号的正交小波变换的分解公式为：其中为尺度系数，为小波系数，h，g是一对正交镜像滤波器组，j为分解层数，N为离散采样点数。从信号滤波的角度看，正交小波分解是将离散信号分别通过一个低通和高通滤波器进行滤波，滤波输出分别对应信号的低频概貌和高频细节，并且每次分解都是对低频部分进行分解，

4、而高频部分则不再继续分解。每次分解的信号长度均减半，相当于在滤波后进行了“二抽一采样”。正交小波变换的重构过程是分解的逆运算，其重构公式为：2.2 小波变换阈值去噪假设N点离散含噪信号为式中为原始信号，为服从的高斯白噪声。对各层小波系数进行阈值处理，传统阈值处理方法有硬阈值法和软阈值法两种，其中硬阈值处理方法是即，比较含噪信号的小波系数与所选定的阈值，小于阈值的点变为0；大于等于阈值的点保持不变。软阈值处理方法是即，比较含噪信号的小波系数与选定的阈值，大于等于阈值的点，收缩为该点值与阈值的差值；小于等于阈值相反数的点，收缩为该点值与阈值的和；阈值小于等于阈值的点变为零。对信号作小波变换后，噪声

5、的小波变换系数主要集中在小尺度上，原始信号的小波变换系数主要集中在大尺度上。小波阈值法消噪就是采用合适的阈值处理方法估计原始信号的小波变换系数并进行重构。2.3 平移不变量小波去噪原理伪吉布斯现象与信号的不连续点的位置有关，更确切地说，和信号的特征（如不连续点）与小波基元素的特征之间的精确对准有关，正交小波变换具有平移不变性，因此可以通过平移含噪信号来改变不连续点的位置，再对平移后的信号进行阈值法去噪处理，然后把去噪后的信号在进行相反的平移，便可以得到原始含噪信号的去噪信号。假如原始含噪信号包含若干个不连续点，它们之间会相互产生干扰，一个不连续点的最佳平移可能是另一个不连续点的最差平移。所以不

6、能采用单一平移，通常采用通过改变平移量，重复平移消噪过程，并把每次平移消噪后的结果求平均的方法消除噪声，即所谓的“平移去噪平均”的平移量小波去噪方法。对于一个信号，我们用表示对信号进行的时域平移，是正整数，即且可逆，令，然后用表示对信号用Donoho的阈值法进行去噪处理，Ave表示平均，则n次循环平移的平移不变量小波去噪方法可以用下式表示： 上式处理过程：首先把原始含噪信号在时域重复平移个单位，再对各次平移后的信号进行小波阈值消噪，然后对各次消噪后的信号进行相同单位的反向平移，最后将各次的处理结果求平均值。3.实验方案本次实验数据来自于“BCI Competiton”实验数据，将数据直接导入m

7、atlab中即可。平移不变量小波去噪具体方法如下：(1)对含噪的EEG信号进行小波分解。根据Mallet最优小波基的选择准则，综合考虑了小波基的对称性、正交性、消失矩及紧支性等数学特性，最终选择Sym8小波基，经对原始含噪信号进行不同分解层数试验比较，最终选择4层分解。因此，本实验采用Sym8小波，对原始含噪信号进行4层小波分解，提取各个尺度上的小波系数。(2)对原始含噪EEG信号进行小波4尺度分解，得到各个尺度上的小波系数，然后选择用软阈值法，分别对各尺度上的高频小波系数进行估计处理。其中，全局阈值，式中为信号长度，由于噪声主要集中在最高分辨率，所以可用小波系数估计，取。(3)结合“平移去噪

8、平均”的平移不变量小波去噪法和阈值量化算法，对含噪EEG信号进行消噪。文中将含噪EEG信号向左循环平移1位，平移4次，然后将各次平移信号消噪后进行相反平移，最后对每次平移后的结果求平均。(4)为了进一步讨论消噪的客观结果，在标准EEG信号中加入高斯白噪声，用Sym8小波对其进行4尺度分解，并用软阈值法和平移不变量小波去噪方法对含噪信号进行消噪。最后用信噪比(SNR)和均方根误差(RMES)来评估消噪结果。信噪比的计算方法是均方根误差的计算方法是实验结果和分析本文所使用的原始EEG信号如下图中的original signal。图(1)，(2)，(3)分别是对噪声信噪比为5dB，15dB，30dB

9、的三种原始加噪信号进行软阈值小波消噪和平移不变量小波消噪的实验仿真图。图(1) 对噪声信噪比为5dB加噪信号消噪图(2) 对噪声信噪比为15dB加噪信号消噪图(3) 对噪声信噪比为30dB加噪信号消噪图(1)(3)只是给出了定性的消噪结果对比分析，为了定量地分析几种消噪方法的结果，本文设计了一种向EEG信号中加入信噪比为5dB、15dB、30dB的高斯白噪声的方法，信号加噪之后再用软阈值法和平移不变量小波去噪分别对其进行消噪处理。最后利用公式计算消噪后信号的信噪比和均方误差，结果如下表：表1 消噪结果的SNR和RMES比较噪声信噪比/dB软阈值法平移不变量小波去噪方法SNRRMESSNRRME

10、S50.895411.99446.05736.6204155.41946.353512.88962.68843015.20552.022921.60680.9681 根据上表所列数据，平移不变量小波去噪方法对加噪的EEG信号消噪结果的信噪比明显改善，且均方误差也明显变小。综合图(1)(3)和表1的分析结果，可以得出平移不变量小波消噪方法对EEG信号的消噪结果优于软阈值法。5.小结本文采用“平移去噪平均”的平移不变量小波去噪的算法对脑电信号进行去噪，经过matlab仿真实验的结果表明，该方法与阈值法相比，能够有效地去除伪吉布斯现象，得到比阈值法去噪更好的视觉效果；用该方法去噪还能减少原信号与去噪

11、信号的均方根，提高信噪比，总的来看，此方法是一种较阈值法去噪更好的一种方法。本文设计的平移不变量小波去噪算法不足之处在于处理高频小波系数采用的软阈值法，在未来学习中应该进一步改进阈值法来处理高频小波系数，以达到更好的去噪效果。附录资料:不需要的可以自行删除c语言典型问题处理方法小结循环问题（1）、数论问题1、求素数 for(i=2;i1，如果它仅有平凡约数1和a，则我们称a为素数(或质数)。整数 1 被称为基数，它既不是质数也不是合数。整数 0 和所有负整数既不是素数，也不是合数。 2、求最大公约数和最小公倍数a、 if(ab) t=a; a=b; b=t; for(i=a;i=1;i-) i

12、f(a%i=0&b%i=0) break; printf(largest common divisor:%dn,i); printf(least common multiple:%dn,(a*b)/is);b、辗转相除法求解 a1=a; b1=b; while(a%b!=0) t=a%b; a=b; b=t; printf(largest common divisor:%dnleast common multiple:%d,b,a1*b1/b);3、求完数 一个数如果恰好等于它的因子之和，这个数就称为“完数”。 例如：6的因子为1、2、3，而6123，因此6是“完数”。for(a=1;a=10

13、00;a+) s=0; for(i=1;i=a) break; if(s=a) printf(%dt,a);注意S=0所放的位置 4、分解质因数 将一个整数写成几个质因数的连乘积，如： 输入36，则程序输出36=2*2*3*3 。解一、看似简单，但要自己完整地写出来还真不容易！竟然还动用了goto语句，正好可以熟悉一下goto语句的用法！main() int a,z,i; clrscr(); scanf(%d,&a);判断下一个数开始有要重新从2开始了。所以用loop语句回到for语句，这是for语句仍从2初始化。从2开始的原则不变，变的是a的值。 loop: for(z=2;z=a;z+)判

14、断是否为质数for(i=2;i=z;i+) if(z%i=0) break;判断是否为a的质因数 if(z=i) if(a%z=0) k+; if(k=1) printf(%d=%d,a1,z);用计数器来解决每行输入不同的问题。 else printf(*%d,z); a/=z; goto loop; 解二：main() int n, k=2, isfirst=1; printf(Input n=); scanf(%d,&n); while(k=n) if(n%k=0) if(isfirst) printf(%d=%d, n, k); isfirst=0; else printf(*%d,k

15、); n/=k; else k+; printf(n);5、从键盘输入两个整数，输出这两个整数的商的小数点后所有1000位整数 for(i=1;i=2;i-) if(fm%i=0&fz%i=0) fz/=i; fm/=i; z=fz/fm; fzx=fz%fm; if(fzx=0) printf(%d%d/%d-%d%d/%d=%dn,z1,fz1,fm1,z2,fz2,fm2,z); else if(z=0) printf(%d%d/%d-%d%d/%d=%d/%dn,z1,fz1,fm1,z2,fz2,fm2,fzx,fm); else printf(%d%d/%d-%d%d/%d=%d%

16、d/%dn,z1,fz1,fm1,z2,fz2,fm2,z,fzx,fm);（2）近似问题1、书P122习题4-6注意千万不要忘记添加#include “math.h”#include math.hmain() float x,j=1,k,s,so; int n; scanf(%f,&x); s=x; so=x+1; for(n=1;fabs(s-so)1e-6;n+) for(k=1;k1e-6) x=(x1+x2)/2; f=x*x*x+4*x*x-10;可以用/*if(f*f10) x2=x; else x1=x; printf(%fn,x);（3）枚举法（4）数列问题二、数组问题（1）

17、排序问题1、从小到大排序main() int a10,i,j,t; for(i=0;i10;i+) scanf(%d,&ai); for(i=1;i10;i+) for(j=0;jaj+1) t=aj+1;aj+1=aj;aj=t; for(i=0;i10;i+) printf(%d ,ai); printf(n);注意排序问题：1、须迅速，熟练，无差错经常插入在程序中间2、现使用最大数下沉冒泡法还可以使用最小数上浮冒泡法3、j控制前面一个数和后面一个数一一比较。由于是最大数下沉，i+1后j仍要从0开始。4、i控制这样的操作一共要做多少次5、注意i j的控制次数2、从大到小排序main()现使

18、用最大数上浮冒泡法还可使用最小数下沉冒泡法 int a10,i,j,t; for(i=0;i10;i+) scanf(%d,&ai); for(i=1;i=i;j-) if(ajaj-1) t=aj-1; aj-1=aj; aj=t; for(i=0;i10;i+) printf(%d ,ai);（2）二维数组三、字符或字符串输入输出问题（1）字符打印1、打印*此类题的溯源为书P122 4.11（1），其他题都是它的拓展 for (i=1;i=n;i+) 一共要输出的行数 for(j=1;j=i;j+) 每行要打印的*数 printf(*); printf(n); a、*解题要点：此类题关键在

19、于找到每行要打印的个数和行数的关系。此题j=i j=n-i+1b、* for(i=1;i=n;i+) 一共要输出的行数 for(j=1;j=n-i;j+) 控制空格数 printf( ); for(k=1;k=i;k+) 每行要打印的*数 printf(*); printf(n); c、 * * *解题要点：在出现空格的时候，在找到每行要打印的*个数和行数的关系后，还应找到空格和行数的关系，分不同的参数进行循环。此题k=i j=n-i j=i-1k=n-i+1d、* * for(i=1;i=n;i+) for(j=1;j=n-i;j+) printf( ); for(k=1;k=2*i-1;k

20、+) printf(*); printf(n); * *e、 * * * for(i=1;i=n-1;i+) for(j=1;j=i;j+) printf( ); for(k=1;k=2*(n-1-i)+1;k+) printf(*); printf(n); for(i=1;i=n;i+) for(j=1;j=n-i;j+) printf( ); for(k=1;k=2*i-1;k+) printf(*); printf(n); * * * * * *2、打印9*9乘法表解题要点：注意寻找行与列的规律。i*ji代表列j代表行for(i=1;i=9;i+) for(j=1;j=9;j+) pri

21、ntf(%-3d ,i*j); 注意输出格式的控制 printf(n); 3、九九乘法表1 2 3 4 5 6 7 8 92 4 6 8 10 12 14 16 183 6 9 12 15 18 21 24 27 9 18 27 36 45 54 63 72 814、杨晖三角形11 11 2 11 3 3 11 4 6 4 11 5 10 10 5 11 6 15 20 15 6 1（2）字符串打印问题for(i=1;i=7;i+) ai1=1; aii=1; for(i=3;i=7;i+) for(j=2;j=i-1;j+) aij=ai-1j-1+ai-1j; gets(a); puts(

22、a); for(i=1;i0;j-) aj=aj-1; a0=t; for(k=0;k=a&ai=z) ai-=32; puts (a);3、逆序输出gets (a); c=strlen(a); for(i=0;i=0;i-) 蓝色部分可以简写为绿色部分coutai-1;4、如输入：ab1 3,;z 输出：ab1注意点：1、= =2、while语句的使用处体会3、全面考虑问题 3,;zgets(a); while(a0= ) for(i=0;ai!=0;i+) ai=ai+1; for(i=0;ai!=0;i+)if(ai= &ai+1!= ) printf(n); else if(ai= &

23、ai+1= ) for(k=i;ak!=0;k+) ak+1=ak+2; i-; elseprintf(%c,ai);5、输入3个字符串，按从小到大排序输出这3个字符串。 使用一个两维数组贮存多个字符串char a8181;注意：如何使用一个两维数组贮存多个字符串 int i,j; for(i=0;i3;i+) gets(ai); for(i=0;i3;i+) puts(ai);注意：1、scanf(%d%s,&n,str) 其中%s为字符串格式2、逐个给字符串赋值的方法见书140页。 不可for（i=0；ai!=0;i+）3、stri=stri-A+10;4、pow函数5、任何进制转为十进制

24、的方法6、输入一个整数n和一个字符串str，计算并输出n进制数str的值。 如输入：7 16则输出：13(16)7=(13)10如输入：16 3A则输出：58(3A)16=(58)10#include stdio.h#include math.hmain() char str81; int n,i,s=0,t; clrscr(); scanf(%d%s,&n,str); for(i=0;stri!=0;i+) if(stri=A) stri=stri-A+10; else stri=stri-0; t=strlen(str); for(i=0;stri!=0;i+) s+=strt-i-1*p

25、ow(n,i); printf(%d,s);编写程序，将一个十进制正整数转换成十六进制数。 注意类比#include main()char a20;int x,i=0,j;clrscr();scanf(%d,&x);while(x) if(x%16=10&x%16=0;j-)printf(%c,aj);printf(n);7、输入一个字符串，将其中的缩写形式展开，并输出展开后的该字符串。所谓展开缩写形式就是将其中由大小写字母或数字构成的形如a-f、U-Z、3-8 的形式展开成为 abcdef 、UVWXYZ 、345678，若出现f-a、A-7、9-5等形式则不予理睬。例如： 输入：qwe24

26、6e-hA-d$-%4-7A-Dz-xp-R4-0输出：qwr246efghA-d$-%4567ABCDz-xp-R4-0main() char a81; int i,c,s,k,t; gets(a); for(i=0;ai!=0;i+) if(ai=-) if(ai-1=A&ai+1=a&ai+1=0&ai+1i;k-)ak+c-2=ak;as-1+c-2+1=0; for(;i=t;i+) ai=ai-1+1; puts(a);补充：循环：求：a+aa+aaa+.的值#includevoid main()int a,n,i=1,sn=0,tn=0;coutinput a and nan;w

27、hile(i=n)tn=tn+a;sn+=tn;a*=10;i+;coutthe answer is snendl;两个乒乓球队进行比赛，各出3人。甲队为A,B,C;已对是X,Y,Z；已经抽签决定比赛名单。有人向队员大厅比赛的名单。A说他不和X比，C说他不和X,Z比。请编程序找出3对赛手的名单。#includevoid main()char i,j,k;for(i=X;i=Z;i+)for(j=X;j=Z;j+)if(i!=j)for (k=X;k=Z;k+)if(i!=k&j!=k)if(i!=X&k!=X&k!=Z)coutA-i B-j C-kendl;枚举口袋中有红，黄，蓝，白，黑5种

28、颜色的球若干。每次从口袋中任意取出3歌，问得到3种不同颜色球的可能取法，输出每种排列的情况。#include#include /在C语言中不用加这句void main()enum colorred ,yellow ,blue,white, black;color pri;int i,j,k,n=0,loop;for(i=red;i=black;i+)for(j=red;j=black;j+)if(i!=j)for (k=red;k=black;k+)if(k!=i)&(k!=j)n+;coutsetw(3)n; /setw是输出格式的限定for(loop=1;loop=3;loop+)swit

29、ch(loop)case 1:pri=color(i);break;case 2:pri=color(j);break;case 3:pri=color(k);break;default:break;switch(pri)case red:coutsetw(8)red;break;case yellow:coutsetw(8)yellow;break;case blue:coutsetw(8)blue;break;case white:coutsetw(8)white;break;case black:coutsetw(8)black;break;default:break;coutendl;

30、couttotal:nendl;数组：不用strcat函数,编写程序，将两个字符串串接起来。#includeusing namespace std;void main()int i=0;char a20,b10;gets(a);gets(b);for(i=0;i10;i+)if(ai=0)break;for(int j=0;j10;j+)ai=bj;i+;if(bj=0)break;coutaendl;输入一行字符，统计其中有多少个单词，单词之间用空格分开。#includeusing namespace std;void main()char string100;int i,num=0,word=0;char c;gets(string);for(i=0;(c=stringi)!=0;i+)if(c= )word=0;else if(word=0)word=1;num+;