如何使用MATLAB和Python實現(xiàn)MFCC特征參數(shù)提取

這篇文章主要介紹如何使用MATLAB和Python實現(xiàn)MFCC特征參數(shù)提取，文中介紹的非常詳細(xì)，具有一定的參考價值，感興趣的小伙伴們一定要看完！

1、MFCC概述

在語音識別（Speech Recognition）和話者識別（Speaker Recognition）方面，最常用到的語音特征就是梅爾倒譜系數(shù)（Mel-scale FrequencyCepstral Coefficients，簡稱MFCC）。根據(jù)人耳聽覺機(jī)理的研究發(fā)現(xiàn)，人耳對不同頻率的聲波有不同的聽覺敏感度。從200Hz到5000Hz的語音信號對語音的清晰度影響較大。兩個響度不等的聲音作用于人耳時，則響度較高的頻率成分的存在會影響到對響度較低的頻率成分的感受，使其變得不易察覺，這種現(xiàn)象稱為掩蔽效應(yīng)。由于頻率較低的聲音在內(nèi)耳蝸基底膜上行波傳遞的距離大于頻率較高的聲音，故一般來說，低音容易掩蔽高音，而高音掩蔽低音較困難。在低頻處的聲音掩蔽的臨界帶寬較高頻要小。所以，人們從低頻到高頻這一段頻帶內(nèi)按臨界帶寬的大小由密到疏安排一組帶通濾波器，對輸入信號進(jìn)行濾波。將每個帶通濾波器輸出的信號能量作為信號的基本特征，對此特征經(jīng)過進(jìn)一步處理后就可以作為語音的輸入特征。由于這種特征不依賴于信號的性質(zhì)，對輸入信號不做任何的假設(shè)和限制，又利用了聽覺模型的研究成果。因此，這種參數(shù)比基于聲道模型的LPCC相比具有更好的魯邦性，更符合人耳的聽覺特性，而且當(dāng)信噪比降低時仍然具有較好的識別性能。

梅爾倒譜系數(shù)是在Mel標(biāo)度頻率域提取出來的倒譜參數(shù)，Mel標(biāo)度描述了人耳頻率的非線性特性，它與頻率的關(guān)系可用下式近似表示：

式中f為頻率，單位為Hz。下圖為Mel頻率與線性頻率的關(guān)系：

2、 MFCC特征參數(shù)提取過程詳解

(1)預(yù)處理

預(yù)處理包括預(yù)加重、分幀、加窗函數(shù)。

預(yù)加重：預(yù)加重的目的是提升高頻部分，使信號的頻譜變得平坦，保持在低頻到高頻的整個頻帶中，能用同樣的信噪比求頻譜。同時，也是為了消除發(fā)生過程中聲帶和嘴唇的效應(yīng)，來補償語音信號受到發(fā)音系統(tǒng)所抑制的高頻部分，也為了突出高頻的共振峰。預(yù)加重處理其實是將語音信號通過一個高通濾波器：

分幀：先將N個采樣點集合成一個觀測單位，稱為幀。通常情況下N的值為256或512，涵蓋的時間約為20~30ms左右。為了避免相鄰兩幀的變化過大，因此會讓兩相鄰幀之間有一段重疊區(qū)域，此重疊區(qū)域包含了M個取樣點，通常M的值約為N的1/2或1/3。通常語音識別所采用語音信號的采樣頻率為8KHz或16KHz，以8KHz來說，若幀長度為256個采樣點，則對應(yīng)的時間長度是256/8000×1000=32ms。

加窗：將每一幀乘以漢明窗，以增加幀左端和右端的連續(xù)性。假設(shè)分幀后的信號為S(n), n=0,1…,N-1, N為幀的大小，那么乘上漢明窗

后 ，W(n)形式如下：

(2)FFT

由于信號在時域上的變換通常很難看出信號的特性，所以通常將它轉(zhuǎn)換為頻域上的能量分布來觀察，不同的能量分布，就能代表不同語音的特性。所以在乘上漢明窗后，每幀還必須再經(jīng)過快速傅里葉變換以得到在頻譜上的能量分布。對分幀加窗后的各幀信號進(jìn)行快速傅里葉變換得到各幀的頻譜。

(3)譜線能量

對語音信號的頻譜取模平方得到語音信號的譜線能量。

(4)計算通過Mel濾波器的能量

將能量譜通過一組Mel尺度的三角形濾波器組，定義一個有M個濾波器的濾波器組（濾波器的個數(shù)和臨界帶的個數(shù)相近），采用的濾波器為三角濾波器，中心頻率為f(m) 。M通常取22-26。各f(m)之間的間隔隨著m值的減小而縮小，隨著m值的增大而增寬，如圖所示：

三角濾波器的頻率響應(yīng)定義為:

對頻譜進(jìn)行平滑化，并消除諧波的作用，突顯原先語音的共振峰。（因此一段語音的音調(diào)或音高，是不會呈現(xiàn)在MFCC 參數(shù)內(nèi)，換句話說，以MFCC 為特征的語音辨識系統(tǒng)，并不會受到輸入語音的音調(diào)不同而有所影響）此外，還可以降低運算量。

計算每個濾波器組輸出的對數(shù)能量為 ：

(5)計算DCT倒譜

經(jīng)離散余弦變換（DCT）得到MFCC系數(shù) :

將上述的對數(shù)能量帶入離散余弦變換，求出L階的Mel參數(shù)。L階指MFCC系數(shù)階數(shù)，通常取12-16。這里M是三角濾波器個數(shù)。

3、MATLAB實現(xiàn)方法

注：在提取MFCC參數(shù)之前需要加載并使用VOICEBOX工具包

Df=5;
fs=8000; 
N=fs/Df;
t=0:1./fs:(N-1)./fs;
x=sin(2*pi*200*t);  
bank=melbankm(24,256,8000,0,0.5,'t');%Mel濾波器的階數(shù)為24，fft變換的長度為256，采樣頻率為8000Hz 
%歸一化mel濾波器組系數(shù) 
bank=full(bank); 
bank=bank/max(bank(:)); 
% DCT系數(shù),12*p
for k=1:12      
n=0:23; 
dctcoef(k,:)=cos((2*n+1)*k*pi/(2*24)); 
end 
%歸一化倒譜提升窗口 
w=1+6*sin(pi*[1:12]./12); 
%w=w/max(w);
%語音信號分幀 
xx=enframe(x,256,80);%對x 256點分為一幀  
%計算每幀的MFCC參數(shù) 
for i=1:size(xx,1) 
y=xx(i,:); 
s=y'.*hamming(256); 
t=abs(fft(s));%fft快速傅立葉變換 
t=t.^2; 
c1=dctcoef*log(bank*t(1:129)); 
c2=c1.*w'; 
end  
plot(c2);title('MFCC');

結(jié)果：

4、Python實現(xiàn)方法

import numpy as np 
from scipy import signal
from scipy.fftpack import dct
import pylab as plt

def enframe(wave_data, nw, inc, winfunc):
  '''將音頻信號轉(zhuǎn)化為幀。
  參數(shù)含義：
  wave_data:原始音頻型號
  nw:每一幀的長度(這里指采樣點的長度，即采樣頻率乘以時間間隔)
  inc:相鄰幀的間隔（同上定義）
  '''
  wlen=len(wave_data) #信號總長度
  if wlen<=nw: #若信號長度小于一個幀的長度，則幀數(shù)定義為1
    nf=1
  else: #否則，計算幀的總長度
    nf=int(np.ceil((1.0*wlen-nw+inc)/inc))
  pad_length=int((nf-1)*inc+nw) #所有幀加起來總的鋪平后的長度
  zeros=np.zeros((pad_length-wlen,)) #不夠的長度使用0填補，類似于FFT中的擴(kuò)充數(shù)組操作
  pad_signal=np.concatenate((wave_data,zeros)) #填補后的信號記為pad_signal
  indices=np.tile(np.arange(0,nw),(nf,1))+np.tile(np.arange(0,nf*inc,inc),(nw,1)).T #相當(dāng)于對所有幀的時間點進(jìn)行抽取，得到nf*nw長度的矩陣
  indices=np.array(indices,dtype=np.int32) #將indices轉(zhuǎn)化為矩陣
  frames=pad_signal[indices] #得到幀信號
  win=np.tile(winfunc,(nf,1)) #window窗函數(shù)，這里默認(rèn)取1
  return frames*win  #返回幀信號矩陣

Df=5
fs=8000
N=fs/Df
t = np.arange(0,(N-1)/fs,1/fs)   
wave_data=np.sin(2*np.pi*200*t)
#預(yù)加重
#b,a = signal.butter(1,1-0.97,'high')
#emphasized_signal = signal.filtfilt(b,a,wave_data)
#歸一化倒譜提升窗口
lifts=[]
for n in range(1,13):
  lift =1 + 6 * np.sin(np.pi * n / 12)
  lifts.append(lift)
#print(lifts)  

#分幀、加窗 
winfunc = signal.hamming(256) 
X=enframe(wave_data, 256, 80, winfunc)  #轉(zhuǎn)置的原因是分幀函數(shù)enframe的輸出矩陣是幀數(shù)*幀長
frameNum =X.shape[0] #返回矩陣行數(shù)18，獲取幀數(shù)
#print(frameNum)
for i in range(frameNum):
  y=X[i,:]
  #fft
  yf = np.abs(np.fft.fft(y)) 
  #print(yf.shape)
  #譜線能量
  yf = yf**2
  #梅爾濾波器系數(shù)
  nfilt = 24
  low_freq_mel = 0
  NFFT=256
  high_freq_mel = (2595 * np.log10(1 + (fs / 2) / 700)) # 把 Hz 變成 Mel
  mel_points = np.linspace(low_freq_mel, high_freq_mel, nfilt + 2) # 將梅爾刻度等間隔
  hz_points = (700 * (10**(mel_points / 2595) - 1)) # 把 Mel 變成 Hz
  bin = np.floor((NFFT + 1) * hz_points / fs)
  fbank = np.zeros((nfilt, int(np.floor(NFFT / 2 + 1))))
  for m in range(1, nfilt + 1):
    f_m_minus = int(bin[m - 1])  # left
    f_m = int(bin[m])       # center
    f_m_plus = int(bin[m + 1])  # right
    for k in range(f_m_minus, f_m):
      fbank[m - 1, k] = (k - bin[m - 1]) / (bin[m] - bin[m - 1])
    for k in range(f_m, f_m_plus):
      fbank[m - 1, k] = (bin[m + 1] - k) / (bin[m + 1] - bin[m])
  filter_banks = np.dot(yf[0:129], fbank.T)
  filter_banks = np.where(filter_banks == 0, np.finfo(float).eps, filter_banks) # 數(shù)值穩(wěn)定性
  filter_banks = 10 * np.log10(filter_banks) # dB 
  filter_banks -= (np.mean(filter_banks, axis=0) + 1e-8)
  #print(filter_banks)
  #DCT系數(shù)
  num_ceps = 12
  c2 = dct(filter_banks, type=2, axis=-1, norm='ortho')[ 1 : (num_ceps + 1)] # Keep 2-13
  c2 *= lifts
print(c2)
plt.plot(c2)
plt.show()

結(jié)果:

由MATLAB和Python繪制出來的波形可以看出二者計算出來的MFCC倒譜系數(shù)的基本走勢相同。

以上是“如何使用MATLAB和Python實現(xiàn)MFCC特征參數(shù)提取”這篇文章的所有內(nèi)容，感謝各位的閱讀！希望分享的內(nèi)容對大家有幫助，更多相關(guān)知識，歡迎關(guān)注億速云行業(yè)資訊頻道！