基于matlab的心音信号处理Word格式.docx
《基于matlab的心音信号处理Word格式.docx》由会员分享,可在线阅读,更多相关《基于matlab的心音信号处理Word格式.docx(20页珍藏版)》请在冰豆网上搜索。
料
1黄文梅,熊佳林,杨勇编著.信号分析与处理——MATALB言及应用.长沙:
国防科技大学出版社,2000
2钱同惠编著.数字信号处理.北京:
机械工业出版社,2004
3姚天任,江太辉编著.数字信号处理.第2版.武汉:
武汉理工大学出版社,2000
指导教师签字
基层教学单位主任签字
说明:
此表一式四份,学生、指导教师、基层教学单位、系部各一份
摘要
心电信号是人体最重要的信号之一,能反映心脏的变时性和变力性,可应用于心血管疾病的诊断和心肌收缩能力的评估;
同时心音信号检测方便、无创、花费极少,可作为心脏疾病检测、预防的有效手段。
因此,研制一种能简易、方便地检测心音信号的数字式心音分析滤波器,对于满足医院和病人的需要,有着极大的社会价值和经济价值。
本课程设计在labview环境下,设计出滤波器编程,可以很好的分析、处理、显示、统计心音信号的信息,充分发挥了微机强大的功能和软件设计的灵活性。
经过运行程序,测试结果显示能够实现从一个包含多种频率成分的心音信号中提取出所需的单一频率心音信号的功能。
关键词:
心电信号
labview
滤波器
1.实验原理
1.1心电产生原理及心电图
我们常说的心电图一般指体表心电图,反映了心脏电兴奋在心脏传导系统中产生和传导的过程。
正常人体的每一个心动周期中,各部分兴奋过程中出现的电变化的方向、途径、次序和时问都有一定的规律,这种生物电变化通过心脏周围的导电组织和体液,反映到身体表面,使身体各部位在每一心动周期中也出现有规律的电变化。
在人体不同部位放置电极,并通过电联线与心电图机的正负极相连,在心电图机上便可以记录到周期变化的心电图。
心电图是通过二次投影形成的。
整体心肌细胞的除极和复极所产生的每一瞬l’日J的除极、复极综合向量轨迹,在立体心脏的三维空『日J内按时问顺序将其顶端相连,便构成立体心向量环。
立体心向量环在额面和横面的投影,形成平面的心向量环;
将平面向量环在导联轴上进行二次投影,就形成相应的心电图。
对于标准的12导联来说,额面心向量环在肢体导联上的投影,形成I、II、Ill、avR、avL、avF导联心电图,而横面心向量环在胸导联轴上的投影便形成了V1~V6导联心电图m。
不同导联记录到的心电图,在波形上有所不同,但基本上都包括一个P波,一个QRS波和一个T波,有时候在T波后还出现一个小u波。
正常的心电图如图卜1所示,下面列出各种波形的定义和生理参数。
图1.1正常心电图
(1)P波:
代表左右心房去极化过程的电位变化。
其波形小而圆钝。
P波在不同导联上的形状可能有差异,在avR导联的P波倒置,在其它的导联中则以向上的正波为主,尤其以II和avF导联最为明显。
历时0.08~O.11s,波幅不超过0.25mv。
(2)QRS波:
代表左、右两心室去极化过程的电位变化。
历时0.06~0.10s。
典型的QRS波包含三个紧密相连的波,第一个向下的波为Q波,其后向上的高而尖的为R波,继R波之后的一个向下的为S波。
但是在不同导联记录的心电图上这三个波不一定都出现,其波形和幅度变化也较大。
(3)T波:
代表心室快速复极化过程中的电位变化。
波形圆钝,历时0.05~0.25s,波形的前肢较长而后肢较短。
T波的方向与QRS波的主波方向一致,在R波为主的导联中,T波的幅度不应低于R波的1/10。
(4)P.R间期:
由P波起点到QRS波起点之问的I’日J期,代表自心房除极开始至心室除极的时『日J。
正常成人为0.12~0.20s。
在幼儿及心率较快的情况下,P.R问期相应缩短;
而经常进行体育锻炼的人,如职业运动员,其P—R问期较长,有的可超过0,20s。
(5)s.T段:
自QRS波终点至T波起点之问的线段,代表心室缓慢复极过程。
正常的ST段多为一等电位线,有时也可有轻微的偏移。
在任一导联,ST段下移一般不应超过0.05mv,ST段上抬在V1,V2导联不超过0.3mv,在V3不超过0.5mv,在V4~V6导联与肢体导联不应超过0.1mv。
(6)Q.T间期:
从QRS波群起点至T波终点的时程,代表心室肌除极和复极全过程所需要的时『自J。
这一时问的长短与心率密切相关。
心率越快,Q.T问期越短;
反之,则问期越长。
心率在60~100次/分时,Q-T怄l期的证常范围为0.32~0,44s。
由于Q-T问期受心率的影响,临床常用校J下的Q.T间期,通常采用Bazett公式计算Q乃=QT/4Rfi。
Q正就是R—Rf.日J期为ls(心律60次/分)时的Q.T问期。
Q正,的正常上限值为0.44s。
(7)U波:
在T波后0.02~0.04s可能出现的低而宽的波,代表心室后继电位。
其方向一般与同导联T波方向一致,幅度较T波低。
u波在肢体导联中不易辨认,一般在胸导联中比较清楚。
1.2首先利用二维plot函数将心电数据生成心电图
①心电数据(170):
表1
心电数据
-2.90
-8.17
-9.14
-7.27
-4.50
-1.78
-0.23
-0.04
0.11
-0.03
-0.15
-0.16
-0.10
0.10
0.70
2.05
5.52
7.26
7.62
5.51
2.29
-0.83
-2.79
10.08
9.41
6.17
2.43
5.12
5.44
1.84
1.57
0.74
-0.65
-0.07
-1.15
-2.94
-0.69
6.23
9.05
9.93
9.50
6.73
3.18
1.86
2.64
2.84
2.33
2.30
1.50
1.11
1.02
1.30
1.10
1.40
3.30
3.09
3.16
3.42
3.63
4.06
4.46
3.92
3.97
2.72
2.52
2.66
4.92
5.46
6.47
5.20
7.89
11.10
1.15
-2.02
0.46
8.07
5.87
15.86
15.98
9.76
1.87
-6.56
-9.62
-6.67
-2.34
0.19
-1.46
-1.58
-1.40
1.19
2.00
4.32
6.14
7.71
7.86
7.15
5.78
3.54
0.31
0.27
12.46
8.68
7.40
7.42
8.79
5.11
2.60
2.28
0.30
-1.27
-2.01
-1.83
-3.79
-1.16
6.05
7.07
6.74
6.70
5.04
2.62
1.91
2.20
1.47
1.64
2.44
3.55
3.87
4.36
4.40
3.82
3.01
1.97
2.09
1.34
0.51
0.29
1.51
2.49
4.51
6.93
6.27
7.95
13.05
4.07
-3.90
1.78
8.66
10.70
3.05
14.35
16.71
13.18
4.59
-3.68
-7.44
-6.36
心电数据存储在maltab\work\ECG.txt文本文件里。
②心电数据生成心电图程序代码:
ECG=load('
ECG.txt'
);
a=length(ECG);
t=[1/a:
1/a:
1]
y=ECG(:
1);
plot(t,y);
title('
含噪心电信号'
xlabel('
时间(ms)'
ylabel('
幅度(dB)'
③心电图:
图1.2心电图
1.3对所生成的心电图进行滤波
一般情况下,我们所要分析的心电数据都存在着不同程度的噪声。
噪声本身就是异常的结构,尤其对于本研究来讲,如果除噪效果不好,很容易就误将噪声识别为奇异结构,这就使得整个研究失去意义。
因此,本研究对心电时间序列的除噪提出了更高的要求。
为了使我们的分析、处理尽量免受噪声的干扰,我们在对心电信号进行分析、处理之前,先对心电数据进行除噪处理。
本模块基于Matlab平台,采用Matlab工具包中提供的小波基函数〔’6〕设计小波滤波器,来提取高频干扰信号;
采用MatIab形态学函数中的开运算和闭运算函数〔’了〕提取基线漂移信号;
并将所提取的高频干扰信号和基线漂移信号作为参考输入,通过编写Matlab代码完成自适应滤波器〔’川的设计,进而完成整个除噪模块的设计。
①小波滤波程序及代码:
%
mallet_wavelet.m
此函数用于研究Mallet算法及滤波器设计
此函数仅用于消噪
a=pi/8;
%角度赋初值
b=pi/8;
%低通重构FIR滤波器h0(n)冲激响应赋值
h0=cos(a)*cos(b);
h1=sin(a)*cos(b);
h2=-sin(a)*sin(b);
h3=cos(a)*sin(b);
low_construct=[h0,h1,h2,h3];
L_fre=4;
%滤波器长度
low_decompose=low_construct(end:
-1:
1);
%确定h0(-n),低通分解滤波器
fori_high=1:
L_fre;
%确定h1(n)=(-1)^n,高通重建滤波器
if(mod(i_high,2)==0);
coefficient=-1;
else
coefficient=1;
end
high_construct(1,i_high)=low_decompose(1,i_high)*coefficient;
high_decompose=high_construct(end:
%高通分解滤波器h1(-n)
L_signal=100;
%信号长度
n=1:
L_signal;
%信号赋值
f=10;
figure
(1);
plot(y);
原信号'
check1=sum(high_decompose);
%h0(n)性质校验
check2=sum(low_decompose);
check3=norm(high_decompose);
check4=norm(low_decompose);
l_fre=conv(y,low_decompose);
%卷积
l_fre_down=dyaddown(l_fre);
%抽取,得低频细节
h_fre=conv(y,high_decompose);
h_fre_down=dyaddown(h_fre);
%信号高频细节
figure
(2);
subplot(2,1,1)
plot(l_fre_down);
小波分解的低频系数'
subplot(2,1,2);
plot(h_fre_down);
小波分解的高频系数'
结果图:
图1.3心电图原图
图1.4小波分解图
②50HZ滤波:
程序代码:
生成一个ideal_lp.m文件
functionhd=ideal_lp(wc,M)
alpha=(M-1)/2;
n=0:
M-1;
m=n-alpha+eps;
%eps为很小的数,避免被0除
hd=sin(wc*m)./(pi*m);
%用Sinc函数产生冲击响应
将ideal_lp.m文件放在matlab\work里
t=[1/(a):
1/(a):
1];
%其中,具有线性相位的FIR低通滤波器由如下函数实现:
%理想低通滤波器
%截止角频率wc,阶数M
%50Hz工频干扰陷波器
%50Hz陷波器:
由一个低通滤波器加上一个高通滤波器组成
%而高通滤波器由一个全通滤波器减去一个低通滤波器构成
M=800;
%滤波器阶数
L=800;
%窗口长度
beta=8;
%衰减系数
Fs=400;
wc1=51/(Fs/2)*pi;
%wc1为高通滤波器截止频率,对应51Hz
wc2=49/(Fs/2)*pi
;
%wc2为低通滤波器截止频率,对应49Hz
h=ideal_lp(pi,M)-ideal_lp(wc1,M)+ideal_lp(wc2,M);
%h为陷波器冲击响应
w=kaiser(L,beta);
b=h.*rot90(w);
%b为50Hz陷波器冲击响应序列
x=filter(b,1,y);
%滤除50Hz工频干扰的心电信号
plot(t,x);
t'
x'
图1.5
50HZ陷波波形
2
结果分析
2.1心电图上的各种波形
一次心动周期就会在新电图上记录出一系列地高低宽窄不同的波形。
包括P波、QRS波群、T波和(无)u波。
了解这些波形及其所代表的意义,是教你怎么看心电图的第二步。
P波,最先出现的一个振幅不高的圆钝波形,它记录的是窦房结激动的右、左心房的激动。
因为窦房结位于右心房,心房的激动先由它开始,所以P波的前半部分记录的是右心房的激动,中间部分记录的是左、右心房的共同激动而后部则代表左心房的激动。
除了aVR导联外,P波基本都是直立的,肢体导联中P波的高度多不超过0.25mV,胸前导联中直立的P波高度不应超过0.15mV。
正常的P波的宽度也不应超过0.11s。
QRS波群,继P波之后出现的一个狭窄但振幅高的波群。
由q波(有或无)、R波和S波组成。
它代表着兴奋从房室结发出先后通过房室束、左右束支和纤细的浦肯野纤维进入心肌细胞,刺激心室的收缩,因此可以将其看作是心室收缩的开始的心电图表现。
Q波,是在出现向上的波之前出现的明确的向下的波形。
如果它很小,宽度不到0.04s,深度不足0.15mV,我们将它记做q波;
若它高且宽,才被称作Q波;
当然有时它是缺无的。
无论有无Q波,第一个出现的向上的高尖的波就是R波;
紧随其后的向下的波就是S波,它也可以根据深度分别命名为S波和s波。
之后出现的向上的波被称作R’(r’)波,向下的波则称作S’(s’)波。
因为波的高低不同,所以可以组合成很多形态,但它也是有限制的,最主要的就是时间限制,通常情况下,正常人的QRS波群的时间0.08s,可以在0.06~0.10s范围内波动。
只要超过这个时限,就应引起注意,特别是超过0.12s便有病理意义了。
T波,上个波群暂停之后出现的波,代表着心室的复极(心室的舒张),以备下一次心室的除极。
观测T波我们要注意它的方向、形态和(高度)深度。
(1)方向,正常情况下,在Ⅰ、Ⅱ导联中T波是直立的;
Ⅲ导联中则可以出现直立、平坦、双向甚至是倒置的T波;
T波在aVR导联中是肯定倒置的,而在aVL和aVF导联中则是和QRS波群的主方向一致的。
胸前导联的T波通常是直立的,当然,V1和V3有时也会出现T波倒置的情况,但它们的深度通常都不会超过0.25mV,当V3导联中出现倒置的T波时,前面两个导连的T波也应该是倒置的,否则就是不正常的表现。
(2)形态,通常T波的波形是圆滑而有个很自然的顶端。
T波一般是不对称的,缓和的上升而略显陡峭地下降至等位线。
(3)高度(深度),各个导联并不完全相同,不过综合看来,在肢体导联中很少超过0.5mV,而在胸前导联中也很少会超过1.0mV。
异常高尖的T波往往出现在心肌梗死的早期或高钾血症。
u波,T波后的一个很微小的波,正常的u波并不是在每一个导联中都显而易见,它究竟代表什么尚无定论。
人体的心电信号极其微弱其幅度一般为几十V至几mV的量级,频率为0.05
Hz~100
Hz。
因此,很容易受到外界的干扰,而且人体作为心电信号的信号源,其内阻是比较大的。
因此在放大部分不但要对心电信号进行放大,还要滤除其他的干扰信号。
图2.1心电图各波段
2.2心电信号噪声分析
心电信号由于受到人体诸多因素的影响,因而有着一般信号所没有的特点:
(1)信号弱,心电信号是体表的电生理信号,一般比较微弱,幅度在10pV~5mV,频率为0.05~100Hz。
例如从母体腹部收取到的胎儿心电信号仅10/zV~50/IV。
(2)噪声强,由于人体自身信号弱,加之人体又是一个复杂的系统,因此信号容易受到噪声干扰。
(3)随机性强,心电信号不仅是随机的,而且是非平稳的。
同时,在心电图检测过程中极易受到各种噪声源的干扰,从而使图像质量变差,使均匀和连续变化的心电数值产生突变,在心电图上形成一些毛刺。
使原本很微弱的信号很难和噪声进行分解。
可能出现的噪声有如下的种类:
(1)工频干扰
工频干扰是由电力系统和人体的分布电容引起的,其频率包括50Hz(MIT-BIH数据库数据工频因为是美国标准,所以是60Hz)的基波及其各次谐波,其幅值成分在ECG峰一峰值的0—50%范围内变化。
(2)引起基线漂移的干扰
心电信号有时候会出现信号基线起伏不平的现象,造成这样的现象有很多原因,主要的有:
①呼吸运动人体呼吸时胸腔内器官和组织会发生一定程度的变化,会对在体表记录到的心电图波形的幅度和形态有所影响,表现为基线随呼吸产生周期性或非周期性漂移,从而导致心电波形的幅度随呼气和吸气而分别上抬和下移。
呼吸运动是引起心电基线漂移的主要原因。
②运动伪迹运动伪迹是由于人体轻微运动造成电极与入体的接触电阻发生变化而引入的一种干扰,它的产生原因仅仅是接触电阻的变化,而不是接触的断续。
这种干扰同样导致信号基线的变化,但不是基线的跃变。
③信号记录和处理中电子设备引起的干扰这种干扰对信号影响很大,严重时可完全淹没心电信号或使得基线剧烈漂移,其中导联开路和放大器的热移是主要因素。
这种干扰往往无法通过心电分析算法来校正。
由于心电波形已经完全畸变,此时对这些数据分析已无太大意义。
所以一般跳过此段数据。
(3)高频噪声
心电信号中的高频噪声主要是肌电噪声。
肌肉收缩会产生mV级的肌电干扰,表现为心电图上不规则的细小波纹,使心电图模糊不清或产生失真。
肌电噪声的特点是频率范围较广,频谱分布非常复杂。
Thakor研究了ECG信号中QRS波、P波和T波以及肌电噪声和运动伪迹的功率谱”
基于小波变换的方法。
近年来一些研究者在QRS波检测中引入了小波变换技术,取得了较好的效果。
小波变换具有良好的时频局域化特性,对时变信号分析有独特的优越性。
通过分析ECG信号功率谱密度的特点以及小波变换和信号频率之问的关系,发现QRS波的能量大多集中在23尺度上,在大于24尺度上则大大减小,而伪迹、基线漂移等能量大多集中于大于25的尺度上。
经过小波变换后,信号的奇异点对应于小波变换的特征值,根据判据规则,找到QRS波的位置,这方面的典型工作是,李翠微等于1995年采用二次样条小波基函数进行21~25尺度的小波变换,采用模极大值算法去检测QRS波m,:
Sahambi提出一个减小高频噪声误判的改进的模极大值算法算法m;
JuanP曲10Martinez等于2004年采用二次样条小波和matlat算法进行QRS波的检测ml;
Tikkanen研究了非线性小波和小波包来消除心电信号的噪声”。
(4)除噪滤波器设计
一般情况下,我们所要分析的心电数据都存在着不同程度的噪声"
噪声本身就是异常的结构,尤其对于本研究来讲,如果除噪效果不好,很容易就误将噪声识别为奇异结构,这就使得整个研究失去意义"
因此,本研究对心电时间序列的除噪提出了更高的要求"
为了使我们的分析!
处理尽量免受噪声的干扰,我们在对心电信号进行分析!
处理之前,先对心电数据进行除噪处理"
。
首先先设计小波滤波器,来提取高频干扰信号,其次采用形态学函数中的开运算和闭运算函数提取基线漂移信号并将所提取的高频干扰信号和基线漂移信号作为参考输入,通过编写代码完成自适应滤波器代码的设计,进而完成整个除噪模块的设计"
常用的除噪方法有小波滤波器!
形态学滤波器和自适应滤波器等等,它们具有良好性能的同时,都还存在着一些缺陷和不足"
如小波滤波器虽然能有效完成信号的滤波处理,但是它在去除低频的基线漂移信号时,必须对信号进行较高尺度的分解与重构,计算量大,实时性差形态滤波器具有计算简单,速度快的特点,对信号中基线漂移的去除,具有近乎完美的表现,但是在滤除高频干扰时会产生信号的高频波动失真而自适应滤波器工作时需要额外采集参考信号,且对参考输入的要求严格,付出的成本高,可靠性降低"
为此,我们在上述几种滤波算法基础上,取长补短,采用基于小波变换与形态学运算的自适应滤波算法1洲,并直接从含噪信号中提取噪声信号作为参考输入"
该算法用小波滤波器提取含噪信号中的高频干扰信号,用形态学滤波器提取低频基线漂移信号,并将这两部分所得到的噪声分量一起作为自适应滤波器的参考输入信号"
通过自适应滤波器的设计和调整,最终获得对信号的满意滤波。
图2.2滤波器设计示意图
3
心得体会
在实验操作中MATLAB程序设计心得:
这次的课程设计让我们对MALTAB有了更加深刻的理解,我们对用MALTAB进行数据分析理解的更加透