comparison toolboxes/MIRtoolbox1.3.2/somtoolbox/som_probability_gmm.m @ 0:e9a9cd732c1e tip

first hg version after svn

+function [pd,Pdm,pmd] = som_probability_gmm(D, sM, K, P)
+
+%SOM_PROBABILITY_GMM Probabilities based on a gaussian mixture model.
+%
+% [pd,Pdm,pmd] = som_probability_gmm(D, sM, K, P)
+% 
+%   [K,P] = som_estimate_gmm(sM,D);
+%   [pd,Pdm,pmd] = som_probability_gmm(D,sM,K,P);
+%   som_show(sM,'color',pmd(:,1),'color',Pdm(:,1))  
+%
+%  Input and output arguments:
+%   D    (matrix) size dlen x dim, the data for which the 
+%        (struct) data struct,     probabilities are calculated
+%   sM   (struct) map struct
+%        (matrix) size munits x dim, the kernel centers
+%   K    (matrix) size munits x dim, kernel width parameters
+%                 computed by SOM_ESTIMATE_GMM
+%   P    (matrix) size 1 x munits, a priori probabilities for each 
+%                 kernel computed by SOM_ESTIMATE_GMM
+%
+%   pd   (vector) size dlen x 1, probability of each data vector in 
+%                 terms of the whole gaussian mixture model
+%   Pdm  (matrix) size munits x dlen, probability of each vector in 
+%                 terms of each kernel
+%   pmd  (matrix) size munits x dlen, probability of each vector to 
+%                 have been generated by each kernel
+%
+% See also SOM_ESTIMATE_GMM.
+
+% Contributed to SOM Toolbox vs2, February 2nd, 2000 by Esa Alhoniemi
+% Copyright (c) by Esa Alhoniemi
+% http://www.cis.hut.fi/projects/somtoolbox/
+
+% ecco 180298 juuso 050100
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+% input arguments
+if isstruct(sM), M = sM.codebook; else M = sM; end
+[c dim] = size(M);
+
+if isstruct(D), D = D.data; end
+dlen = size(D,1);
+
+% reserve space for output variables
+pd = zeros(dlen,1); 
+if nargout>=2, Pdm = zeros(c,dlen); end
+if nargout==3, pmd = zeros(c,dlen); end
+
+% the parameters of each kernel
+cCoeff = cell(c,1);
+cCoinv = cell(c,1);
+for m=1:c, 
+  co = diag(K(m,:));
+  cCoinv{m} = inv(co);
+  cCoeff{m} = 1 / ((2*pi)^(dim/2)*det(co)^.5);
+end
+
+% go through the vectors one by one
+for i=1:dlen, 
+
+  x = D(i,:);
+  
+  % compute p(x|m)
+  pxm = zeros(c,1); 
+  for m = 1:c,
+    dx     = M(m,:) - x;
+    pxm(m) = cCoeff{m} * exp(-.5 * dx * cCoinv{m} * dx');
+    %pxm(m) = normal(dx, zeros(1,dim), diag(K(m,:)));  
+  end
+  pxm(isnan(pxm(:))) = 0;
+  
+  % p(x|m)  
+  if nargin>=2, Pdm(:,i) = pxm; end
+  
+  % P(x) = P(x|M) = sum( P(m) * p(x|m) )
+  pd(i) = P*pxm; 
+  
+  % p(m|x) = p(x|m) * P(m) / P(x)
+  if nargout==3, pmd(:,i) = (P' .* pxm) / pd(i); end
+  
+end
+
+
+return; 
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% subfunction normal
+%
+% computes probability of x when mean and covariance matrix
+% of a distribution are known
+
+function result = normal(x, mu, co)
+
+[l dim] = size(x);
+coinv   = inv(co);
+coeff   = 1 / ((2*pi)^(dim/2)*det(co)^.5);
+diff   = x - mu;
+result = coeff * exp(-.5 * diff * coinv * diff');