声音模型的识别与建立-精品资料.doc

1、2016 长沙理工大学第二届研究生数学建模竞赛参赛承诺书我们仔细阅读了长沙理工大学研究生数学建模竞赛的竞赛规则.我们完全明白，在竞赛开始后参赛队员不能以任何方式 （包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。我们知道，抄袭别人的成果是违反竞赛规则的 , 如果引用别人的成果或其他公开的资料（包括网上查到的资料） ，必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。我们授权长沙理工大学研究生数学建模竞赛组委会，可将我们的论文以任何形式进

2、行公开展示（包括进行网上公示，在书籍、期刊和其他媒体进行正式或非正式发表等）。我们参赛选择的题号是（从组委会提供的试题中选择一项填写）：我们的参赛报名号为（如果组委会设置报名号的话）：所属学校（请填写完整的全名）：参赛队员 ( 打印并签名 ) ： 1.2.3.指导教师或指导教师组负责人( 打印并签名 ) ：日期：年月日评阅编号（由组委会评阅前进行编号）：2016 长沙理工大学研究生数学建模竞赛编号专用页评阅编号（由组委会评阅前进行编号）：评阅记录（可供评阅时使用）：评阅人评分备注声音识别模型的建立与评价摘 要本文通过使用 MATLAB软件对声音的时域和频域特征进行了提取，研究特征向量提取方法及

3、 SVM核函数和参数选取对识别结果的影响， 分析特征提取算法的优缺点， 以及不同核函数以及惩罚参数对识别性能的影响。 通过使用支持向量机法建模， 基本达到了区分正常声音与非正常声音的目的。最后提出用低通滤波滤除白噪声。关键词： 特征向量，支持向量机，核函数，低通滤波，白噪声从物理上讲，声音是由物体振动产生的一种波，并通过空气作用于人的耳鼓，使人们能够感知。声音的具有四种性质：1）音高：振动发出的声波有不同的频率，称为“音高”；2）强弱：声音的强弱是由振幅决定的，振幅是代表物体振动强度的特定单位，一般用分贝 (dB)来表示。 3）长短：一般把声音的发展过程分为四个阶段，分别是触发、衰减、保持和消

4、失。这四个阶段称为“包络” ，包络的发生时间，也就是一个声音的长短。4) 音质：音质好的声音听起来悦耳，相反则让人不适。问题一利用 matlab 中的 sound 函数，播放出声音信号，试听并比较正常和非正常开门声音的差别，利用 plot 函数绘制出具体的声音波形图，总结差别在哪些方面？试听：我们使用 sound 函数播放声音样本。 在仔细听了正常开门及非正常开门的声音后，发现了他们之间的差别：正常开门声音很短促，即声音的长短度短，且其强弱度相对较低。相反非正常开门声音持续时间长，强弱度高。画图：我们选取了三组正常开门及三组非正常开门的声音，使用 plot 函数画出声音波形图，如下图 1 所示

5、：由图 1 可以看出正常开门的声音波形比较疏松，所以音调就低。非正常开门的声音波形比较密集，所以声调就高。 我们还可以看出前者声波比较集中， 而后者则比较散， 即跨度大。这一点很好理解， 正常用钥匙开门所需时间肯定比盗贼撬锁所需时间短， 所以就造成了这种现象。程序源代码见附录一。图 1 正常与非正常开门声音波形图问题二利用合适的时域或（和）频域特征表达个声音信号，建立特征向量，写出提取特征向量的具体方法和程序代码。首先，我们对两种样本求其均值。由于正常开门的第一组数据有人声干扰，故舍去。然后用 plot 函数画出声音波形图，如下图 2 所示：根据所分析的参数类型， 语音信号分析可以分成时域分析

6、和变换域( 频域、倒谱域 ) 分析。其中时域分析方法是最简单、最直观的方法，因为它直接对语音信号的时域波形进行分析。接下来我们先进行时域分析。图 2 正常开门声与非正常开门声均值短时平均能量：定义n 时刻某语音信号的短时平均能量En 为：nEn( )(n) 2x mmmm n ( N 1)()(n) 2x mm（ 1）式中， N 为窗长，可见短时能量为一帧样点值的平方和。一般我们认为声音在 10-30ms 之内是稳定的，取帧长也在 10-30ms 之内，而帧移通常取 5-15ms 之间，所以取 N=55、95、125、165。如下图 3、4 为正常开门及非正常开门 N 取不同值时短时能量函数随

7、帧数的变化曲线，其中横坐标为帧数。由图 3、图 4 可以看出， N=55,N=95时的曲线不够平滑，而 N=165的曲线又过于平滑，故选取 N=125时的曲线。图 3 正常开门声音的短时能量曲线图 4 非正常开门声音的短时能量曲线通过观察短时能量曲线，可以看出正常开门时的能量比较集中且数值小，非正常开门时的能量比较分散且数值大。容易想到，这些现象与前面问题一的结论是相吻合的。短时平均过零率：短时平均过零率是指每帧内信号通过零值的次数。对有时间横轴的连续声音信号， 可以观察到声音的时域波形通过横轴的情况。 在离散时间声音信号情况下， 如果相邻的采样具有不同的代数符号就称为发生了过零，因此可以计算

8、过零的次数。Zn| sgnx( m) sgnx( m1)| ( n m)m（2）= |sgnx( n) sgnx( n 1) *( n)上式为短时平均过零率的公式，其中，sgn 为符号函数，即sgn x( n)1,x( n)01,x( n)0（ 3）短时平均过零率曲线如下图5、图 6 所示：图 5 正常开门声音的短时平均过零率图 6 非正常开门声音的短时平均过零率由图 5、图 6 可以看出前者比后者的过零率要高。短时自相关函数：自相关函数用于衡量信号自身时间波形的相似性。对于声音来说，采用短时分析方法，可以定义短时自相关函数为：Rnk)x mn mx m k) (n k m（4）( ) ()

9、()m短时自相关函数曲线如下图7、图 8 所示：图 7 正常开门声音的自相关函数曲线图 8 非正常开门声音的自相关函数频域特征：频域分析主要是对声音波形进行傅里叶变换，如下图9 所示：图9声音波形的傅里叶变换由图9 可以看出正常开门的对数幅度要比非正常开门的对数幅度小的多。程序源代码见附录二。问题三建立声音识别模型（二分类模型） ，利用模型区分正常和非正常声音，评价模型的好坏。问题三利用支持向量机构建二分类器对正常与非正常敲门进行分类预测。 支持向量机是 Cortes 和 Vapnik 于 1995 年提出的，与传统分类器比较，它在解决小样本、非线性及高维模式识别中表现出许多特有的优势， 并能

10、够推广应用到函数拟合等其他机器学习问题中。 支持向量机方法是建立在统计学习理论的 VC 维理论和结构风险最小原理基础上的， 根据有限的样本信息在模型的复杂性和学习能力之间寻求最佳折衷，以期获得最好的推广能力。图 10 支持向量机分类器结构图采集样本共有 80 组音频信号，从中选择 60 组数据做为训练数据训练网络， 20 组数据作为测试数据测试网络分类能力。下面是准备数据步骤：Step1:根据给定的数据，选定训练集和测试集；Step2:为训练集与测试集选定标签集；Step3:利用训练集进行训练分类器得到model；Step4:根据 model，对测试集进行测试集得到accuracy rate；

11、训练集： trainset();分别取正常开门与非正常开门数据的一半作为训练集；测试集： testset();分别取正常开门与非正常开门数据的另一半作为测试集；标签集： label();取 bedroom的数据为正类标签为1；forse 的数据为负类标签为 -1.下面为 svm二分类器的分类结果：图 11 基于支持向量机的分类器上图将展示了训练与测试的结果， 结果等于 1 表示正常声音，结果为 -1 为非正常开门声，通过学习训练得到上图结果。程序源代码见附件三。通过仿真，验证了支持向量机在解决分类问题中的优势，今后可以将这方法用于更深层次的语音识别中。问题四试利用特征选择或变换，对特征向量进行

12、优化，并利用参数优化技术优化模型的参数，使识别模型的准确率提高。问题四其中支持向量机的核函数主要有：线性核函数，多项式核函数，径向基核函数，两层感知器核函数。不同的核函数直接影响SVM的分类能力，使用不同的核函数SVM分类的效果不同，所需要的时间也不同，下图是使用短时能量值和均值为特征参数，选取RBF为核函数的分类效果图：图 12 SVMRBF分类图图 13 分类模型错误率此模型相比于问题三中建立的 SVM分类模型分类正确率明显提高， 达到了 85.65%，可见合适的核函数选取的重要意义。程序源代码见附录四。问题五白噪声，是指功率谱密度在整个频域内均匀分布的噪声，即其功率谱密度为pnN 0 2

13、(5)式中， N 0 为常数，。由于 1 和 (t ) 为一对傅里叶变换对，所以白噪声的自相关函数为RnN 0(6)2由式 (1) 、式(2) 可知，白噪声的自相关函数仅在0 时才不为零 , 而对于其他任意的，白噪声的自相关函数都为零，即在任意两个不同时刻上的随机变量都是不相关的。对于一个线性时不变系统 ( 如滤波器 ) ，该系统是由它的冲激响应 h(t) 或等效地由它的频率响应 H ( f ) 表征，这里的 h(t ) 和 H ( f ) 是一对傅里叶变换对。 如果令 x(t ) 为系统的输人信号，y(t) 是输出信号。系统的输出可以表示成如下形式：y thx td(7)如果 x(t) 是平

14、稳随机过程 X (t ) 的样本函数，那么 y(t) 就是随机过程 Y(t ) 的样本函数。可以求出输出的自相关函数是RhhRdd(8)yyxx由于我们知道自相关函数和功率谱密度函数是一对傅里叶变换对，所以可以得到输出过程的功率密度谱，即为相关函数的傅里叶变换：j 2fxx fH f2(9)yy fRyy ed由此可以看出，输出信号的功率谱密度就是输入信号的功率谱密度乘以系统的频率响应的模的平方。当输入随机过程是白噪声时， 输出随机过程的自相关特性和功率密度谱将完全由系统的频率响应所决定。我们设计了一个低通滤波器，对加了白噪声的声音信号进行滤波，即可去除白噪声。低通滤波 (Low-pass f

15、ilter) 是一种过滤方式，规则为低频信号能正常通过，而超过设定临界值的高频信号则被阻隔、 减弱。但是阻隔、 减弱的幅度则会依据不同的频率以及不同的滤波程序（目的）而改变。在数字图像处理领域，从频域看，低通滤波可以对图像进行平滑去噪处理。滤波器将这部分不需要的频率滤除，就实现了去噪的效果。 滤波前与滤波后的声音波形对比如下图所示：图 14 滤波前后信号波形对比实验表明低通滤波具有良好的去噪能力，该方法对添加的噪声处理后，与原始信号效果几乎拟合。证明低通滤波对信号去噪是可行的。程序源代码见附录五附 录附录一：clc;clearfor i=5:7j=2*(i-4)-1;a=正 num2str(i

16、);load(a)subplot(3,2,j);plot(y)title(正 num2str(i)xlabel(Time);ylabel(Amplitude);a=非 num2str(i);load(a)subplot(3,2,j+1)plot(y)title(非 num2str(i)xlabel(Time);ylabel(Amplitude);endsaveas(gcf,q1.jpg);附录二：均值function fei=feipj()fei=zeros(55125,1);for i=1:40a=非 num2str(i);load(a)fei=fei+y;endfei=fei/40;短时平

17、均能量function energy=nengliang(y)s=fra(55,55,y);s2=s.2;energy=sum(s2,2)%;%一帧内各样点的能量求一帧能量subplot(2,2,1)%定义画图数量和布局plot(energy)%画 N=55时的语音能量图xlabel(帧数 )%横坐标ylabel(短时能量 E)%纵坐标title(N=55)%曲线标识s=fra(95,95,y) ;s2=s.2;energy=sum(s2,2);subplot(2,2,2)plot(energy)%画 N=95时的语音能量图xlabel(帧数 )ylabel(短时能量 E)title(N=95

18、)s=fra(125,125,y);s2=s.2;energy=sum(s2,2);subplot(2,2,3)plot(energy)%画 N=125时的语音能量图xlabel(帧数 )ylabel(短时能量 E)title(N=125)s=fra(165,165,y);s2=s.2;energy=sum(s2,2);subplot(2,2,4)plot(energy)%画 N=165时的语音能量图xlabel(帧数 )ylabel(短时能量 E)title(N=165)function f=fra(len,inc,x)fh=fix(size(x,1)-len)/inc)+1);f=zero

19、s(fh,len);i=1;n=1;while i=fhj=1;while j=lenf(i,j)=x(n);j=j+1;n=n+1;endn=n-len+inc;i=i+1;end短时平均过零率function f=zcro(x)f=zeros(size(x,1),1);%生成全零矩阵for i=1:size(x,1)z=x(i,:);%提取一行数据for j=1:(length(z)-1);if z(j)*z(j+1)0;f(i)=f(i)+1;endendendfunction zcr=guoling(x)s=fra(125,125,x);%分帧，帧移 110zcr=zcro(s);%求

20、过零率figure(1);subplot(2,1,1)plot(x);title(非);xlabel(Time);ylabel(Amplitude);subplot(2,1,2)plot(zcr);xlabel(帧数 );ylabel(过零次数 );title(原始信号的过零率 );短时自相关函数function A1=zixiangguan(s1)N=125;A=;%加 N=125的矩形窗for k=1:125;sum=0;for m=1:N-k+1;sum=sum+s1(m)*s1(m+k-1);%计算自相关endA(k)=sum;endfor k=1:125A1(k)=A(k)/A(1)

21、;%归一化 A(k);endplot(A1);title(非正常开门的自相关函数)xlabel(延时 k)ylabel(R(k)傅里叶变换clc;clearFs=11025;zheng=zhengpj();fei=feipj();T=1/Fs;L=55125;t=(0:L-1)*T;NFFT=2nextpow2(L);subplot(2,1,1)Y = fft(zheng,NFFT)/L;f = Fs/2*linspace(0,1,NFFT/2+1);y1=2*abs(Y(1:NFFT/2+1);plot(f,2*abs(Y(1:NFFT/2+1)title(正常开门的傅里叶变换 )subpl

22、ot(2,1,2)Y = fft(fei,NFFT)/L;f = Fs/2*linspace(0,1,NFFT/2+1);y2=2*abs(Y(1:NFFT/2+1);plot(f,2*abs(Y(1:NFFT/2+1)title(非正常开门的傅里叶变换)saveas(gcf,fuliye.jpg)附录三：clc;clearzhn=zeros(40,55125);zhn(1,:)=zhengpj();for i=2:40name=正,num2str(i);load(name)zhn(i,:)=y;endfei=zeros(40,55125);for i=1:40name=非,num2str(i

23、);load(name)fei(i,:)=y;endlabel=ones(1,641),-ones(1,641);trainset=zhn(1:20,:);fei(1:20,:);testset=zhn(21:40,:);fei(21:40,:);mtrain,ntrain = size(trainset);mtest,ntest = size(testset);test_dataset = trainset;testset;dataset_scale,ps = mapminmax(test_dataset,0,1);dataset_scale = dataset_scale;trainset

24、 = dataset_scale(1:mtrain,:);testset = dataset_scale( (mtrain+1):(mtrain+mtest),: ); train=zeros(2,55125);test=zeros(2,55125);for i=1:20train(1,:)=trainset(i,:)+train(1,:);train(2,:)=trainset(i+20,:)+train(2,:);endfor i=1:20test(1,:)=testset(i,:)+test(1,:);test(2,:)=testset(i+20,:)+test(2,:);endtrai

25、n1=zeros(2,1282);test1=zeros(2,1282);train1(1,:)=train(1,1:43:55125);train1(2,:)=train(2,1:43:55125);test1(1,:)=test(1,1:43:55125);test1(2,:)=test(2,1:43:55125);train1=train1;test1=test1;p=cvpartition(label,Holdout,0.2);model = svmtrain(train1(p.training,:),label(p.training), showplot,true); C = svm

26、classify(model,train1(p.test,:),showplot,true);err_rate = sum(label(p.test)= C)/p.TestSizesaveas(gcf,wentisan2.jpg);附录四：clear;clc;closeall;disp(-disp(输入样本矩阵采用径向基内积函数的支持向量机的应用( 每行表示一个测试数据 ):);-);%x=randn(1,20)-3randn(1,20)+3;randn(1,40)loadpj1.matloadpj2.matdisp(样本所属类别:);a1=pj1(1:100:55125);a2=pj2(1:100:55125);a3=a1a2;y=-ones