聲音模型的識別與建立-精品資料.doc

1、2016 長沙理工大學(xué)第二屆研究生數(shù)學(xué)建模競賽參賽承諾書我們仔細閱讀了長沙理工大學(xué)研究生數(shù)學(xué)建模競賽的競賽規(guī)則.我們完全明白，在競賽開始后參賽隊員不能以任何方式 （包括電話、電子郵件、網(wǎng)上咨詢等）與隊外的任何人（包括指導(dǎo)教師）研究、討論與賽題有關(guān)的問題。我們知道，抄襲別人的成果是違反競賽規(guī)則的 , 如果引用別人的成果或其他公開的資料（包括網(wǎng)上查到的資料） ，必須按照規(guī)定的參考文獻的表述方式在正文引用處和參考文獻中明確列出。我們鄭重承諾，嚴格遵守競賽規(guī)則，以保證競賽的公正、公平性。如有違反競賽規(guī)則的行為，我們將受到嚴肅處理。我們授權(quán)長沙理工大學(xué)研究生數(shù)學(xué)建模競賽組委會，可將我們的論文以任何形式進

2、行公開展示（包括進行網(wǎng)上公示，在書籍、期刊和其他媒體進行正式或非正式發(fā)表等）。我們參賽選擇的題號是（從組委會提供的試題中選擇一項填寫）：我們的參賽報名號為（如果組委會設(shè)置報名號的話）：所屬學(xué)校（請?zhí)顚懲暾娜簠①愱爢T ( 打印并簽名 ) ： 1.2.3.指導(dǎo)教師或指導(dǎo)教師組負責人( 打印并簽名 ) ：日期：年月日評閱編號（由組委會評閱前進行編號）：2016 長沙理工大學(xué)研究生數(shù)學(xué)建模競賽編號專用頁評閱編號（由組委會評閱前進行編號）：評閱記錄（可供評閱時使用）：評閱人評分備注聲音識別模型的建立與評價摘 要本文通過使用 MATLAB軟件對聲音的時域和頻域特征進行了提取，研究特征向量提取方法及

3、 SVM核函數(shù)和參數(shù)選取對識別結(jié)果的影響， 分析特征提取算法的優(yōu)缺點， 以及不同核函數(shù)以及懲罰參數(shù)對識別性能的影響。 通過使用支持向量機法建模， 基本達到了區(qū)分正常聲音與非正常聲音的目的。最后提出用低通濾波濾除白噪聲。關(guān)鍵詞： 特征向量，支持向量機，核函數(shù)，低通濾波，白噪聲從物理上講，聲音是由物體振動產(chǎn)生的一種波，并通過空氣作用于人的耳鼓，使人們能夠感知。聲音的具有四種性質(zhì)：1）音高：振動發(fā)出的聲波有不同的頻率，稱為“音高”；2）強弱：聲音的強弱是由振幅決定的，振幅是代表物體振動強度的特定單位，一般用分貝 (dB)來表示。 3）長短：一般把聲音的發(fā)展過程分為四個階段，分別是觸發(fā)、衰減、保持和消

4、失。這四個階段稱為“包絡(luò)” ，包絡(luò)的發(fā)生時間，也就是一個聲音的長短。4) 音質(zhì)：音質(zhì)好的聲音聽起來悅耳，相反則讓人不適。問題一利用 matlab 中的 sound 函數(shù)，播放出聲音信號，試聽并比較正常和非正常開門聲音的差別，利用 plot 函數(shù)繪制出具體的聲音波形圖，總結(jié)差別在哪些方面？試聽：我們使用 sound 函數(shù)播放聲音樣本。 在仔細聽了正常開門及非正常開門的聲音后，發(fā)現(xiàn)了他們之間的差別：正常開門聲音很短促，即聲音的長短度短，且其強弱度相對較低。相反非正常開門聲音持續(xù)時間長，強弱度高。畫圖：我們選取了三組正常開門及三組非正常開門的聲音，使用 plot 函數(shù)畫出聲音波形圖，如下圖 1 所示

5、：由圖 1 可以看出正常開門的聲音波形比較疏松，所以音調(diào)就低。非正常開門的聲音波形比較密集，所以聲調(diào)就高。 我們還可以看出前者聲波比較集中， 而后者則比較散， 即跨度大。這一點很好理解， 正常用鑰匙開門所需時間肯定比盜賊撬鎖所需時間短， 所以就造成了這種現(xiàn)象。程序源代碼見附錄一。圖 1 正常與非正常開門聲音波形圖問題二利用合適的時域或（和）頻域特征表達個聲音信號，建立特征向量，寫出提取特征向量的具體方法和程序代碼。首先，我們對兩種樣本求其均值。由于正常開門的第一組數(shù)據(jù)有人聲干擾，故舍去。然后用 plot 函數(shù)畫出聲音波形圖，如下圖 2 所示：根據(jù)所分析的參數(shù)類型， 語音信號分析可以分成時域分析

6、和變換域( 頻域、倒譜域 ) 分析。其中時域分析方法是最簡單、最直觀的方法，因為它直接對語音信號的時域波形進行分析。接下來我們先進行時域分析。圖 2 正常開門聲與非正常開門聲均值短時平均能量：定義n 時刻某語音信號的短時平均能量En 為：nEn( )(n) 2x mmmm n ( N 1)()(n) 2x mm（ 1）式中， N 為窗長，可見短時能量為一幀樣點值的平方和。一般我們認為聲音在 10-30ms 之內(nèi)是穩(wěn)定的，取幀長也在 10-30ms 之內(nèi)，而幀移通常取 5-15ms 之間，所以取 N=55、95、125、165。如下圖 3、4 為正常開門及非正常開門 N 取不同值時短時能量函數(shù)隨

7、幀數(shù)的變化曲線，其中橫坐標為幀數(shù)。由圖 3、圖 4 可以看出， N=55,N=95時的曲線不夠平滑，而 N=165的曲線又過于平滑，故選取 N=125時的曲線。圖 3 正常開門聲音的短時能量曲線圖 4 非正常開門聲音的短時能量曲線通過觀察短時能量曲線，可以看出正常開門時的能量比較集中且數(shù)值小，非正常開門時的能量比較分散且數(shù)值大。容易想到，這些現(xiàn)象與前面問題一的結(jié)論是相吻合的。短時平均過零率：短時平均過零率是指每幀內(nèi)信號通過零值的次數(shù)。對有時間橫軸的連續(xù)聲音信號， 可以觀察到聲音的時域波形通過橫軸的情況。 在離散時間聲音信號情況下， 如果相鄰的采樣具有不同的代數(shù)符號就稱為發(fā)生了過零，因此可以計算

8、過零的次數(shù)。Zn| sgnx( m) sgnx( m1)| ( n m)m（2）= |sgnx( n) sgnx( n 1) *( n)上式為短時平均過零率的公式，其中，sgn 為符號函數(shù)，即sgn x( n)1,x( n)01,x( n)0（ 3）短時平均過零率曲線如下圖5、圖 6 所示：圖 5 正常開門聲音的短時平均過零率圖 6 非正常開門聲音的短時平均過零率由圖 5、圖 6 可以看出前者比后者的過零率要高。短時自相關(guān)函數(shù)：自相關(guān)函數(shù)用于衡量信號自身時間波形的相似性。對于聲音來說，采用短時分析方法，可以定義短時自相關(guān)函數(shù)為：Rnk)x mn mx m k) (n k m（4）( ) ()

9、()m短時自相關(guān)函數(shù)曲線如下圖7、圖 8 所示：圖 7 正常開門聲音的自相關(guān)函數(shù)曲線圖 8 非正常開門聲音的自相關(guān)函數(shù)頻域特征：頻域分析主要是對聲音波形進行傅里葉變換，如下圖9 所示：圖9聲音波形的傅里葉變換由圖9 可以看出正常開門的對數(shù)幅度要比非正常開門的對數(shù)幅度小的多。程序源代碼見附錄二。問題三建立聲音識別模型（二分類模型） ，利用模型區(qū)分正常和非正常聲音，評價模型的好壞。問題三利用支持向量機構(gòu)建二分類器對正常與非正常敲門進行分類預(yù)測。 支持向量機是 Cortes 和 Vapnik 于 1995 年提出的，與傳統(tǒng)分類器比較，它在解決小樣本、非線性及高維模式識別中表現(xiàn)出許多特有的優(yōu)勢， 并能

10、夠推廣應(yīng)用到函數(shù)擬合等其他機器學(xué)習(xí)問題中。 支持向量機方法是建立在統(tǒng)計學(xué)習(xí)理論的 VC 維理論和結(jié)構(gòu)風(fēng)險最小原理基礎(chǔ)上的， 根據(jù)有限的樣本信息在模型的復(fù)雜性和學(xué)習(xí)能力之間尋求最佳折衷，以期獲得最好的推廣能力。圖 10 支持向量機分類器結(jié)構(gòu)圖采集樣本共有 80 組音頻信號，從中選擇 60 組數(shù)據(jù)做為訓(xùn)練數(shù)據(jù)訓(xùn)練網(wǎng)絡(luò)， 20 組數(shù)據(jù)作為測試數(shù)據(jù)測試網(wǎng)絡(luò)分類能力。下面是準備數(shù)據(jù)步驟：Step1:根據(jù)給定的數(shù)據(jù)，選定訓(xùn)練集和測試集；Step2:為訓(xùn)練集與測試集選定標簽集；Step3:利用訓(xùn)練集進行訓(xùn)練分類器得到model；Step4:根據(jù) model，對測試集進行測試集得到accuracy rate；

11、訓(xùn)練集： trainset();分別取正常開門與非正常開門數(shù)據(jù)的一半作為訓(xùn)練集；測試集： testset();分別取正常開門與非正常開門數(shù)據(jù)的另一半作為測試集；標簽集： label();取 bedroom的數(shù)據(jù)為正類標簽為1；forse 的數(shù)據(jù)為負類標簽為 -1.下面為 svm二分類器的分類結(jié)果：圖 11 基于支持向量機的分類器上圖將展示了訓(xùn)練與測試的結(jié)果， 結(jié)果等于 1 表示正常聲音，結(jié)果為 -1 為非正常開門聲，通過學(xué)習(xí)訓(xùn)練得到上圖結(jié)果。程序源代碼見附件三。通過仿真，驗證了支持向量機在解決分類問題中的優(yōu)勢，今后可以將這方法用于更深層次的語音識別中。問題四試利用特征選擇或變換，對特征向量進行

12、優(yōu)化，并利用參數(shù)優(yōu)化技術(shù)優(yōu)化模型的參數(shù)，使識別模型的準確率提高。問題四其中支持向量機的核函數(shù)主要有：線性核函數(shù)，多項式核函數(shù)，徑向基核函數(shù)，兩層感知器核函數(shù)。不同的核函數(shù)直接影響SVM的分類能力，使用不同的核函數(shù)SVM分類的效果不同，所需要的時間也不同，下圖是使用短時能量值和均值為特征參數(shù)，選取RBF為核函數(shù)的分類效果圖：圖 12 SVMRBF分類圖圖 13 分類模型錯誤率此模型相比于問題三中建立的 SVM分類模型分類正確率明顯提高， 達到了 85.65%，可見合適的核函數(shù)選取的重要意義。程序源代碼見附錄四。問題五白噪聲，是指功率譜密度在整個頻域內(nèi)均勻分布的噪聲，即其功率譜密度為pnN 0 2

13、(5)式中， N 0 為常數(shù)，。由于 1 和 (t ) 為一對傅里葉變換對，所以白噪聲的自相關(guān)函數(shù)為RnN 0(6)2由式 (1) 、式(2) 可知，白噪聲的自相關(guān)函數(shù)僅在0 時才不為零 , 而對于其他任意的，白噪聲的自相關(guān)函數(shù)都為零，即在任意兩個不同時刻上的隨機變量都是不相關(guān)的。對于一個線性時不變系統(tǒng) ( 如濾波器 ) ，該系統(tǒng)是由它的沖激響應(yīng) h(t) 或等效地由它的頻率響應(yīng) H ( f ) 表征，這里的 h(t ) 和 H ( f ) 是一對傅里葉變換對。 如果令 x(t ) 為系統(tǒng)的輸人信號，y(t) 是輸出信號。系統(tǒng)的輸出可以表示成如下形式：y thx td(7)如果 x(t) 是平

14、穩(wěn)隨機過程 X (t ) 的樣本函數(shù)，那么 y(t) 就是隨機過程 Y(t ) 的樣本函數(shù)?？梢郧蟪鲚敵龅淖韵嚓P(guān)函數(shù)是RhhRdd(8)yyxx由于我們知道自相關(guān)函數(shù)和功率譜密度函數(shù)是一對傅里葉變換對，所以可以得到輸出過程的功率密度譜，即為相關(guān)函數(shù)的傅里葉變換：j 2fxx fH f2(9)yy fRyy ed由此可以看出，輸出信號的功率譜密度就是輸入信號的功率譜密度乘以系統(tǒng)的頻率響應(yīng)的模的平方。當輸入隨機過程是白噪聲時， 輸出隨機過程的自相關(guān)特性和功率密度譜將完全由系統(tǒng)的頻率響應(yīng)所決定。我們設(shè)計了一個低通濾波器，對加了白噪聲的聲音信號進行濾波，即可去除白噪聲。低通濾波 (Low-pass f

15、ilter) 是一種過濾方式，規(guī)則為低頻信號能正常通過，而超過設(shè)定臨界值的高頻信號則被阻隔、 減弱。但是阻隔、 減弱的幅度則會依據(jù)不同的頻率以及不同的濾波程序（目的）而改變。在數(shù)字圖像處理領(lǐng)域，從頻域看，低通濾波可以對圖像進行平滑去噪處理。濾波器將這部分不需要的頻率濾除，就實現(xiàn)了去噪的效果。 濾波前與濾波后的聲音波形對比如下圖所示：圖 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)%;%一幀內(nèi)各樣點的能量求一幀能量subplot(2,2,1)%定義畫圖數(shù)量和布局plot(energy)%畫 N=55時的語音能量圖xlabel(幀數(shù) )%橫坐標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(幀數(shù) )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(幀數(shù) )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(幀數(shù) )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,:);%提取一行數(shù)據(jù)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(幀數(shù) );ylabel(過零次數(shù) );title(原始信號的過零率 );短時自相關(guān)函數(shù)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);%計算自相關(guān)endA(k)=sum;endfor k=1:125A1(k)=A(k)/A(1)

21、;%歸一化 A(k);endplot(A1);title(非正常開門的自相關(guān)函數(shù))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(輸入樣本矩陣采用徑向基內(nèi)積函數(shù)的支持向量機的應(yīng)用( 每行表示一個測試數(shù)據(jù) ):);-);%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