?

GMM-HMM語音識別

轉自@http://www.cnblogs.com/tornadomeet

?　　為了對GMM-HMM在語音識別上的應用有個宏觀認識，花了些時間讀了下HTK（用htk完成簡單的孤立詞識別）的部分源碼，對該算法總算有了點大概認識，達到了預期我想要的。不得不說，網絡上關于語音識別的通俗易懂教程太少，都是各種公式滿天飛，很少有說具體細節的，當然了，那需要有實戰經驗才行。下面總結以下幾點，對其有個宏觀印象即可（以孤立詞識別為例）。

　　一、每個單詞的讀音都對應一個HMM模型，大家都知道HMM模型中有個狀態集S，那么每個狀態用什么來表示呢，數字？向量？矩陣？其實這個狀態集中的狀態沒有具體的數學要求，只是一個名稱而已，你可以用’1’, ’2’, ‘3’…表示，也可以用’a’, ‘b’, ’c ’表示。另外每個HMM模型中到底該用多少個狀態，是通過先驗知識人為設定的。

　　二、HMM的每一個狀態都對應有一個觀察值，這個觀察值可以是一個實數，也可以是個向量，且每個狀態對應的觀察值的維度應該相同。假設現在有一個單詞的音頻文件，首先需要將其進行采樣得到數字信息（A/D轉換），然后分幀進行MFCC特征提取，假設每一幀音頻對應的MFCC特征長度為39，則每個音頻文件就轉換成了N個MFCC向量（不同音頻文件對應的N可能不同），這就成了一個序列，而在訓練HMM模型的參數時（比如用Baum-Welch算法），每次輸入到HMM中的數據要求就是一個觀測值序列。這時，每個狀態對應的觀測值為39維的向量，因為向量中元素的取值是連續的，需要用多維密度函數來模擬，通常情況下用的是多維高斯函數。在GMM-HMM體系中，這個擬合函數是用K個多維高斯混合得到的。假設知道了每個狀態對應的K個多維高斯的所有參數，則該GMM生成該狀態上某一個觀察向量（一幀音頻的MFCC系數）的概率就可以求出來了。

　　三、對每個單詞建立一個HMM模型，需要用到該單詞的訓練樣本，這些訓練樣本是提前標注好的，即每個樣本對應一段音頻，該音頻只包含這個單詞的讀音。當有了該單詞的多個訓練樣本后，就用這些樣本結合Baum-Welch算法和EM算法來訓練出GMM-HMM的所有參數，這些參數包括初始狀態的概率向量，狀態之間的轉移矩陣，每個狀態對應的觀察矩陣（這里對應的是GMM，即每個狀態對應的K個高斯的權值，每個高斯的均值向量和方差矩陣）。

　　四、在識別階段，輸入一段音頻，如果該音頻含有多個單詞，則可以手動先將其分割開（考慮的是最簡單的方法），然后提取每個單詞的音頻MFCC特征序列，將該序列輸入到每個HMM模型（已提前訓練好的）中，采用前向算法求出每個HMM模型生成該序列的概率，最后取最大概率對應的那個模型，而那個模型所表示的單詞就是我們識別的結果。

　　五、在建立聲學模型時，可以用Deep Learning的方法來代替GMM-HMM中的GMM，因為GMM模擬任意函數的功能取決于混合高斯函數的個數，所以具有一定的局限性，屬于淺層模型。而Deep Network可以模擬任意的函數，因而表達能力更強。注意，這里用來代替GMM的Deep Nets模型要求是產生式模型，比如DBN，DBM等，因為在訓練HMM-DL網絡時，需要用到HMM的某個狀態產生一個樣本的概率。

　　六、GMM-HMM在具體實現起來還是相當復雜的。

　　七、一般涉及到時間序列時才會使用HMM，比如這里音頻中的語音識別，視頻中的行為識別等。如果我們用GMM-HMM對靜態的圖片分類，因為這里沒涉及到時間信息，所以HMM的狀態數可設為1，那么此時的GMM-HMM算法就退化成GMM算法了。

?

　　MFCC:

　　MFCC的matlab實現教程可參考：張智星老師的網頁教程mfcc.?最基本的12維特征。

function mfcc=frame2mfcc(frame, fs, filterNum, mfccNum, plotOpt) % frame2mfcc: Frame to MFCC conversion. % Usage: mfcc=frame2mfcc(frame, fs, filterNum, mfccNum, plotOpt) % % For example: % waveFile='what_movies_have_you_seen_recently.wav'; % [y, fs, nbits]=wavReadInt(waveFile); % startIndex=12000; % frameSize=512; % frame=y(startIndex:startIndex+frameSize-1); % frame2mfcc(frame, fs, 20, 12, 1);% Roger Jang 20060417if nargin<1, selfdemo; return; end if nargin<2, fs=16000; end if nargin<3, filterNum=20; end if nargin<4, mfccNum=12; end if nargin<5, plotOpt=0; endframeSize=length(frame); % ====== Preemphasis should be done at wave level %a=0.95; %frame2 = filter([1, -a], 1, frame); frame2=frame; % ====== Hamming windowing frame3=frame2.*hamming(frameSize); % ====== FFT [fftMag, fftPhase, fftFreq, fftPowerDb]=fftOneSide(frame3, fs); % ====== Triangular band-pass filter bank triFilterBankPrm=getTriFilterBankPrm(fs, filterNum); % Get parameters for triangular band-pass filter bank % Triangular bandpass filter. for i=1:filterNumtbfCoef(i)=dot(fftPowerDb, trimf(fftFreq, triFilterBankPrm(:,i)));%得到filterNum個濾波系數 end % ====== DCT mfcc=zeros(mfccNum, 1); %DCT變換的前后個數也沒有變 for i=1:mfccNumcoef = cos((pi/filterNum)*i*((1:filterNum)-0.5))'; %mfcc中的前mfccNum個系數mfcc(i) = sum(coef.*tbfCoef');%直接按照DCT公式 end % ====== Log energy %logEnergy=10*log10(sum(frame.*frame)); %mfcc=[logEnergy; mfcc];if plotOptsubplot(2,1,1);plot(frame, '.-');set(gca, 'xlim', [-inf inf]);title('Input frame');subplot(2,1,2);plot(mfcc, '.-');set(gca, 'xlim', [-inf inf]);title('MFCC vector'); end% ====== trimf.m (from fuzzy toolbox) function y = trimf(x, prm) %由頻率的橫坐標算出三角形內的縱坐標,0~1 a = prm(1); b = prm(2); c = prm(3); y = zeros(size(x)); % Left and right shoulders (y = 0) index = find(x <= a | c <= x); y(index) = zeros(size(index)); %只考慮三角波內的量 % Left slope if (a ~= b)index = find(a < x & x < b);y(index) = (x(index)-a)/(b-a); end % right slope if (b ~= c)index = find(b < x & x < c);y(index) = (c-x(index))/(c-b); end % Center (y = 1) index = find(x == b); y(index) = ones(size(index));% ====== Self demo function selfdemo waveFile='what_movies_have_you_seen_recently.wav'; [y, fs, nbits]=wavReadInt(waveFile); startIndex=12000; frameSize=512; frame=y(startIndex:startIndex+frameSize-1); feval(mfilename, frame, fs, 20, 12, 1);

?

　　ZCR:

　　過0檢測，用于判斷每一幀中過零點的數量情況，最簡單的版本可參考：zeros cross rate.?

waveFile='csNthu.wav'; frameSize=256; overlap=0; [y, fs, nbits]=wavread(waveFile); frameMat=enframe(y, frameSize, overlap); frameNum=size(frameMat, 2); for i=1:frameNumframeMat(:,i)=frameMat(:,i)-mean(frameMat(:,i)); % mean justification end zcr=sum(frameMat(1:end-1, :).*frameMat(2:end, :)<0); sampleTime=(1:length(y))/fs; frameTime=((0:frameNum-1)*(frameSize-overlap)+0.5*frameSize)/fs; subplot(2,1,1); plot(sampleTime, y); ylabel('Amplitude'); title(waveFile); subplot(2,1,2); plot(frameTime, zcr, '.-'); xlabel('Time (sec)'); ylabel('Count'); title('ZCR');

?

　　EPD:

　　端點檢測，檢測聲音的起始點和終止點，可參考：EPD in Time Domain,在時域中的最簡單檢測方法。

waveFile='sunday.wav'; [wave, fs, nbits] = wavread(waveFile); frameSize = 256; overlap = 128;wave=wave-mean(wave); % zero-mean substraction frameMat=buffer2(wave, frameSize, overlap); % frame blocking,每一列代表一幀 frameNum=size(frameMat, 2); % no. of frames volume=frame2volume(frameMat); % volume,求每一幀的能量，絕對值或者平方和,volume為行向量 volumeTh1=max(volume)*0.1; % volume threshold 1 volumeTh2=median(volume)*0.1; % volume threshold 2 volumeTh3=min(volume)*10; % volume threshold 3 volumeTh4=volume(1)*5; % volume threshold 4 index1 = find(volume>volumeTh1); %找出volume大于閾值的那些幀序號 index2 = find(volume>volumeTh2); index3 = find(volume>volumeTh3); index4 = find(volume>volumeTh4); %frame2sampleIndex()為從幀序號找到樣本點的序號(即每一個采樣點的序號) %endPointX長度為2,包含了起點和終點的樣本點序號 endPoint1=frame2sampleIndex([index1(1), index1(end)], frameSize, overlap); endPoint2=frame2sampleIndex([index2(1), index2(end)], frameSize, overlap); endPoint3=frame2sampleIndex([index3(1), index3(end)], frameSize, overlap); endPoint4=frame2sampleIndex([index4(1), index4(end)], frameSize, overlap);subplot(2,1,1); time=(1:length(wave))/fs; plot(time, wave); ylabel('Amplitude'); title('Waveform'); axis([-inf inf -1 1]); line(time(endPoint1( 1))*[1 1], [-1, 1], 'color', 'm');%標起點終點線 line(time(endPoint2( 1))*[1 1], [-1, 1], 'color', 'g'); line(time(endPoint3( 1))*[1 1], [-1, 1], 'color', 'k'); line(time(endPoint4( 1))*[1 1], [-1, 1], 'color', 'r'); line(time(endPoint1(end))*[1 1], [-1, 1], 'color', 'm'); line(time(endPoint2(end))*[1 1], [-1, 1], 'color', 'g'); line(time(endPoint3(end))*[1 1], [-1, 1], 'color', 'k'); line(time(endPoint4(end))*[1 1], [-1, 1], 'color', 'r'); legend('Waveform', 'Boundaries by threshold 1', 'Boundaries by threshold 2', 'Boundaries by threshold 3', 'Boundaries by threshold 4');subplot(2,1,2); frameTime=frame2sampleIndex(1:frameNum, frameSize, overlap); plot(frameTime, volume, '.-'); ylabel('Sum of Abs.'); title('Volume'); axis tight; line([min(frameTime), max(frameTime)], volumeTh1*[1 1], 'color', 'm'); line([min(frameTime), max(frameTime)], volumeTh2*[1 1], 'color', 'g'); line([min(frameTime), max(frameTime)], volumeTh3*[1 1], 'color', 'k'); line([min(frameTime), max(frameTime)], volumeTh4*[1 1], 'color', 'r'); legend('Volume', 'Threshold 1', 'Threshold 2', 'Threshold 3', 'Threshold 4');

?

　 ?GMM:?

　 ?GMM用在擬合數據分布上，本質上是先假設樣本的概率分布為GMM，然后用多個樣本去學習這些GMM的參數。GMM建模在語音中可用于某個單詞的發音，某個人的音色等。其訓練過程可參考:speaker recognition.

function [M, V, W, logProb] = gmmTrain(data, gaussianNum, dispOpt) % gmmTrain: Parameter training for gaussian mixture model (GMM) % Usage: function [M, V, W, logProb] = gmm(data, gaussianNum, dispOpt) % data: dim x dataNum matrix where each column is a data point % gaussianNum: No. of Gaussians or initial centers % dispOpt: Option for displaying info during training % M: dim x meanNum matrix where each column is a mean vector % V: 1 x gaussianNum vector where each element is a variance for a Gaussian % W: 1 x gaussianNum vector where each element is a weighting factor for a Gaussian% Roger Jang 20000610if nargin==0, selfdemo; return; end if nargin<3, dispOpt=0; endmaxLoopCount = 50; % Max. iteration minImprove = 1e-6; % Min. improvement minVariance = 1e-6; % Min. variance logProb = zeros(maxLoopCount, 1); % Array for objective function [dim, dataNum] = size(data);% Set initial parameters % Set initial M %M = data(1+floor(rand(gaussianNum,1)*dataNum),:); % Randomly select several data points as the centers if length(gaussianNum)==1,% Using vqKmeans to find initial centersfprintf('Start KMEANS to find the initial mu...\n'); % M = vqKmeansMex(data, gaussianNum, 0);M = vqKmeans(data, gaussianNum, 0); %利用聚類的方法求均值,聚成gaussianNum類 % M = vqLBG(data, gaussianNum, 0);fprintf('Start GMM training...\n');if any(any(~isfinite(M))); keyboard; end else% gaussianNum is in fact the initial centersM = gaussianNum;gaussianNum = size(M, 2); end % Set initial V as the distance to the nearest center if gaussianNum==1V=1; elsedistance=pairwiseSqrDist(M);%pairwiseSqrDist是dll%distance=pairwiseSqrDist2(M);distance(1:(gaussianNum+1):gaussianNum^2)=inf; % Diagonal elements are inf[V, index]=min(distance); % Initial variance for each Gaussian end % Set initial W W = ones(1, gaussianNum)/gaussianNum; % Weight for each Gaussian,初始化時是均分權值if dispOpt & dim==2, displayGmm(M, V, data); end for i = 1:maxLoopCount %開始迭代訓練參數,EM算法% Expectation step:% P(i,j) is the probability of data(:,j) to the i-th Gaussian% Prob為每個樣本在GMM下的概率[prob, P]=gmmEval(data, M, V, W);logProb(i)=sum(log(prob)); %所有樣本的聯合概率if dispOptfprintf('i = %d, log prob. = %f\n',i-1, logProb(i));endPW = diag(W)*P;BETA=PW./(ones(gaussianNum,1)*sum(PW)); % BETA(i,j) is beta_i(x_j)sumBETA=sum(BETA,2);% Maximization step: eqns (2.96) to (2.98) from Bishop p.67:M = (data*BETA')./(ones(dim,1)*sumBETA');DISTSQ = pairwiseSqrDist(M, data); % Distance of M to data%DISTSQ = pairwiseSqrDist2(M, data); % Distance of M to dataV = max((sum(BETA.*DISTSQ, 2)./sumBETA)/dim, minVariance); % (2.97)W = (1/dataNum)*sumBETA; % (2.98)if dispOpt & dim==2, displayGmm(M, V, data); endif i>1, if logProb(i)-logProb(i-1)<minImprove, break; end; end end [prob, P]=gmmEval(data, M, V, W); logProb(i)=sum(log(prob)); fprintf('Iteration count = %d, log prob. = %f\n',i, logProb(i)); logProb(i+1:maxLoopCount) = [];% ====== Self Demo ====== function selfdemo %[data, gaussianNum] = dcdata(2); data = rand(1000,2); gaussianNum = 8; data=data'; plotOpt=1; [M, V, W, lp] = feval(mfilename, data, gaussianNum, plotOpt);pointNum = 40; x = linspace(min(data(1,:)), max(data(1,:)), pointNum); y = linspace(min(data(2,:)), max(data(2,:)), pointNum); [xx, yy] = meshgrid(x, y); data = [xx(:) yy(:)]'; z = gmmEval(data, M, V, W); zz = reshape(z, pointNum, pointNum); figure; mesh(xx, yy, zz); axis tight; box on; rotate3d on figure; contour(xx, yy, zz, 30); axis image% ====== Other subfunctions ====== function displayGmm(M, V, data) % Display function for EM algorithm figureH=findobj(0, 'tag', mfilename); if isempty(figureH)figureH=figure;set(figureH, 'tag', mfilename);colordef blackplot(data(1,:), data(2,:),'.r'); axis imagetheta=linspace(-pi, pi, 21);x=cos(theta); y=sin(theta);sigma=sqrt(V);for i=1:length(sigma)circleH(i)=line(x*sigma(i)+M(1,i), y*sigma(i)+M(2,i), 'color', 'y');endset(circleH, 'tag', 'circleH', 'erasemode', 'xor'); elsecircleH=findobj(figureH, 'tag', 'circleH');theta=linspace(-pi, pi, 21);x=cos(theta); y=sin(theta);sigma=sqrt(V);for i=1:length(sigma)set(circleH(i), 'xdata', x*sigma(i)+M(1,i), 'ydata', y*sigma(i)+M(2,i));enddrawnow end

?

?　?Speaker identification:

　 ?給N個人的語音資料，用GMM可以訓練這N個人的聲音模型，然后給定一段語音，判斷該語音與這N個人中哪個最相似。方法是求出該語音在N個GMM模型下的概率，選出概率最大的那個。可參考:speaker recognition.

function [recogRate, confusionMatrix, speakerData]=speakerIdentify(speakerData, speakerGmm, useIntGmm) % speakerIdentify: speaker identification using GMM parameters % Usage: [recogRate, confusionMatrix, speakerData]=speakerIdentify(speakerData, speakerGmm, useIntGmm) % speakerData: structure array generated by speakerDataRead.m % speakerGmm: speakerGmm(i).gmmPrm is the GMM parameters for speaker i. % useIntGmm: use fixed-point GMM% Roger Jang, 20070517, 20080726if nargin<3, useIntGmm=0; end% ====== Speaker identification using GMM parameters speakerNum=length(speakerData); for i=1:speakerNum % fprintf('%d/%d: Recognizing wave files by %s\n', i, speakerNum, speakerData(i).name);for j=1:length(speakerData(i).sentence) % fprintf('\tSentece %d...\n', j);frameNum=size(speakerData(i).sentence(j).fea, 2);logProb=zeros(speakerNum, frameNum); %logProb(i,m)表示第i個人第j個句子中第m幀在GMM模型下的log概率%找出一個句子，看它屬于哪個speakerfor k=1:speakerNum, % fprintf('\t\tSpeaker %d...\n', k);% logProb(k, :)=gmmEval(speakerData(i).sentence(j).fea, speakerGmm(k).gmmPrm);if ~useIntGmm% logProb(k, :)=gmmEvalMex(speakerData(i).sentence(j).fea, gmm(k).mean, gmm(k).covariance, gmm(k).weight);logProb(k, :)=gmmEval(speakerData(i).sentence(j).fea, speakerGmm(k).gmmPrm);else% logProb(k, :)=gmmEvalIntMex(speakerData(i).sentence(j).fea, gmm(k).mean, gmm(k).covariance, gmm(k).weight);logProb(k, :)=gmmEvalIntMex(speakerData(i).sentence(j).fea, speakerGmm(i).gmmPrm);endendcumLogProb=sum(logProb, 2);[maxProb, index]=max(cumLogProb);speakerData(i).sentence(j).predictedSpeaker=index; %找出身份speakerData(i).sentence(j).logProb=logProb;end end% ====== Compute confusion matrix and recognition rate confusionMatrix=zeros(speakerNum); for i=1:speakerNum,predictedSpeaker=[speakerData(i).sentence.predictedSpeaker];[index, count]=elementCount(predictedSpeaker);confusionMatrix(i, index)=count; end recogRate=sum(diag(confusionMatrix))/sum(sum(confusionMatrix));

?

　　GMM-HMM:

　　訓練階段：給出HMM的k個狀態，每個狀態下的觀察樣本的生成可以用一個概率分布來擬合，這里是采用GMM擬合的。其實，可以把GMM-HMM整體看成是一個生成模型。給定該模型的5個初始參數(結合隨機和訓練樣本獲得)，啟動EM算法的E步：獲得訓練樣本分布，即計算訓練樣本在各個狀態下的概率。M步：用這些訓練樣本重新評估那5個參數。

　　測試階段：(以孤立詞識別為例)給定每個詞發音的frame矩陣，取出某一個GMM-HMM模型，算出該發音每一幀數據在取出的GMM-HMM模型各個state下的概率，結合模型的轉移概率和初始概率，獲得對應的clique tree，可用圖模型的方法inference出生成該語音的概率。比較多個GMM-HMM模型，取最大概率的模型對應的詞。

?

　??參考資料：

? ? ?機器學習&數據挖掘筆記_13（用htk完成簡單的孤立詞識別）

? ? ?http://htk.eng.cam.ac.uk/extensions/

? ? ?張智星老師的網頁教程mfcc.

?

?

總結

以上是生活随笔為你收集整理的GMM-HMM语音识别的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。