hmm.m
资源名称:speech.rar [点击查看]
上传用户:ay_070428
上传日期:2014-12-04
资源大小:11427k
文件大小:3k
源码类别:
语音合成与识别
开发平台:
Matlab
- % function [Mu,Cov,P,Pi,LL]=hmm(X,T,K,cyc,tol);
- %
- % Gaussian Observation Hidden Markov Model
- %
- % X - N x p data matrix
- % T - length of each sequence (N must evenly divide by T, default T=N)
- % K - number of states (default 2)
- % cyc - maximum number of cycles of Baum-Welch (default 100)
- % tol - termination tolerance (prop change in likelihood) (default 0.0001)
- %
- % Mu - mean vectors
- % Cov - output covariance matrix (full, tied across states)
- % P - state transition matrix
- % Pi - priors
- % LL - log likelihood curve
- %
- % Iterates until a proportional change < tol in the log likelihood
- % or cyc steps of Baum-Welch
- function [Mu,Cov,P,Pi,LL]=hmm(X,T,K,cyc,tol)
- p=length(X(1,:));
- N=length(X(:,1));
- if nargin<5 tol=0.0001; end;
- if nargin<4 cyc=100; end;
- if nargin<3 K=2; end;
- if nargin<2 T=N; end;
- if (rem(N,T)~=0)
- disp('Error: Data matrix length must be multiple of sequence length T');
- return;
- end;
- N=N/T;
- Cov=diag(diag(cov(X)));
- Mu=randn(K,p)*sqrtm(Cov)+ones(K,1)*mean(X);
- Pi=rand(1,K);
- Pi=Pi/sum(Pi);
- P=rand(K);
- P=rdiv(P,rsum(P));
- LL=[];
- lik=0;
- alpha=zeros(T,K);
- beta=zeros(T,K);
- gamma=zeros(T,K);
- B=zeros(T,K);
- k1=(2*pi)^(-p/2);
- for cycle=1:cyc
- %%%% FORWARD-BACKWARD
- Gamma=[];
- Gammasum=zeros(1,K);
- Scale=zeros(T,1);
- Xi=zeros(T-1,K*K);
- for n=1:N
- iCov=inv(Cov);
- k2=k1/sqrt(det(Cov));
- for i=1:T
- for l=1:K
- d=Mu(l,:)-X((n-1)*T+i,:);
- B(i,l)=k2*exp(-0.5*d*iCov*d');
- end;
- end;
- scale=zeros(T,1);
- alpha(1,:)=Pi.*B(1,:);
- scale(1)=sum(alpha(1,:));
- alpha(1,:)=alpha(1,:)/scale(1);
- for i=2:T
- alpha(i,:)=(alpha(i-1,:)*P).*B(i,:);
- scale(i)=sum(alpha(i,:));
- alpha(i,:)=alpha(i,:)/scale(i);
- end;
- beta(T,:)=ones(1,K)/scale(T);
- for i=T-1:-1:1
- beta(i,:)=(beta(i+1,:).*B(i+1,:))*(P')/scale(i);
- end;
- gamma=(alpha.*beta);
- gamma=rdiv(gamma,rsum(gamma));
- gammasum=sum(gamma);
- xi=zeros(T-1,K*K);
- for i=1:T-1
- t=P.*( alpha(i,:)' * (beta(i+1,:).*B(i+1,:)));
- xi(i,:)=t(:)'/sum(t(:));
- end;
- Scale=Scale+log(scale);
- Gamma=[Gamma; gamma];
- Gammasum=Gammasum+gammasum;
- Xi=Xi+xi;
- end;
- %%%% M STEP
- % outputs
- Mu=zeros(K,p);
- Mu=Gamma'*X;
- Mu=rdiv(Mu,Gammasum');
- % transition matrix
- sxi=rsum(Xi')';
- sxi=reshape(sxi,K,K);
- P=rdiv(sxi,rsum(sxi));
- % priors
- Pi=zeros(1,K);
- for i=1:N
- Pi=Pi+Gamma((i-1)*T+1,:);
- end
- Pi=Pi/N;
- % covariance
- Cov=zeros(p,p);
- for l=1:K
- d=(X-ones(T*N,1)*Mu(l,:));
- Cov=Cov+rprod(d,Gamma(:,l))'*d;
- end;
- Cov=Cov/(sum(Gammasum));
- oldlik=lik;
- lik=sum(Scale);
- LL=[LL lik];
- fprintf('cycle %i log likelihood = %f ',cycle,lik);
- if (cycle<=2)
- likbase=lik;
- elseif (lik<oldlik)
- fprintf('violation');
- elseif ((lik-likbase)<(1 + tol)*(oldlik-likbase)|~finite(lik))
- fprintf('n');
- break;
- end;
- fprintf('n');
- end