|
wolffd@0
|
1 function [pd,Pdm,pmd] = som_probability_gmm(D, sM, K, P)
|
|
wolffd@0
|
2
|
|
wolffd@0
|
3 %SOM_PROBABILITY_GMM Probabilities based on a gaussian mixture model.
|
|
wolffd@0
|
4 %
|
|
wolffd@0
|
5 % [pd,Pdm,pmd] = som_probability_gmm(D, sM, K, P)
|
|
wolffd@0
|
6 %
|
|
wolffd@0
|
7 % [K,P] = som_estimate_gmm(sM,D);
|
|
wolffd@0
|
8 % [pd,Pdm,pmd] = som_probability_gmm(D,sM,K,P);
|
|
wolffd@0
|
9 % som_show(sM,'color',pmd(:,1),'color',Pdm(:,1))
|
|
wolffd@0
|
10 %
|
|
wolffd@0
|
11 % Input and output arguments:
|
|
wolffd@0
|
12 % D (matrix) size dlen x dim, the data for which the
|
|
wolffd@0
|
13 % (struct) data struct, probabilities are calculated
|
|
wolffd@0
|
14 % sM (struct) map struct
|
|
wolffd@0
|
15 % (matrix) size munits x dim, the kernel centers
|
|
wolffd@0
|
16 % K (matrix) size munits x dim, kernel width parameters
|
|
wolffd@0
|
17 % computed by SOM_ESTIMATE_GMM
|
|
wolffd@0
|
18 % P (matrix) size 1 x munits, a priori probabilities for each
|
|
wolffd@0
|
19 % kernel computed by SOM_ESTIMATE_GMM
|
|
wolffd@0
|
20 %
|
|
wolffd@0
|
21 % pd (vector) size dlen x 1, probability of each data vector in
|
|
wolffd@0
|
22 % terms of the whole gaussian mixture model
|
|
wolffd@0
|
23 % Pdm (matrix) size munits x dlen, probability of each vector in
|
|
wolffd@0
|
24 % terms of each kernel
|
|
wolffd@0
|
25 % pmd (matrix) size munits x dlen, probability of each vector to
|
|
wolffd@0
|
26 % have been generated by each kernel
|
|
wolffd@0
|
27 %
|
|
wolffd@0
|
28 % See also SOM_ESTIMATE_GMM.
|
|
wolffd@0
|
29
|
|
wolffd@0
|
30 % Contributed to SOM Toolbox vs2, February 2nd, 2000 by Esa Alhoniemi
|
|
wolffd@0
|
31 % Copyright (c) by Esa Alhoniemi
|
|
wolffd@0
|
32 % http://www.cis.hut.fi/projects/somtoolbox/
|
|
wolffd@0
|
33
|
|
wolffd@0
|
34 % ecco 180298 juuso 050100
|
|
wolffd@0
|
35
|
|
wolffd@0
|
36 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
|
wolffd@0
|
37
|
|
wolffd@0
|
38 % input arguments
|
|
wolffd@0
|
39 if isstruct(sM), M = sM.codebook; else M = sM; end
|
|
wolffd@0
|
40 [c dim] = size(M);
|
|
wolffd@0
|
41
|
|
wolffd@0
|
42 if isstruct(D), D = D.data; end
|
|
wolffd@0
|
43 dlen = size(D,1);
|
|
wolffd@0
|
44
|
|
wolffd@0
|
45 % reserve space for output variables
|
|
wolffd@0
|
46 pd = zeros(dlen,1);
|
|
wolffd@0
|
47 if nargout>=2, Pdm = zeros(c,dlen); end
|
|
wolffd@0
|
48 if nargout==3, pmd = zeros(c,dlen); end
|
|
wolffd@0
|
49
|
|
wolffd@0
|
50 % the parameters of each kernel
|
|
wolffd@0
|
51 cCoeff = cell(c,1);
|
|
wolffd@0
|
52 cCoinv = cell(c,1);
|
|
wolffd@0
|
53 for m=1:c,
|
|
wolffd@0
|
54 co = diag(K(m,:));
|
|
wolffd@0
|
55 cCoinv{m} = inv(co);
|
|
wolffd@0
|
56 cCoeff{m} = 1 / ((2*pi)^(dim/2)*det(co)^.5);
|
|
wolffd@0
|
57 end
|
|
wolffd@0
|
58
|
|
wolffd@0
|
59 % go through the vectors one by one
|
|
wolffd@0
|
60 for i=1:dlen,
|
|
wolffd@0
|
61
|
|
wolffd@0
|
62 x = D(i,:);
|
|
wolffd@0
|
63
|
|
wolffd@0
|
64 % compute p(x|m)
|
|
wolffd@0
|
65 pxm = zeros(c,1);
|
|
wolffd@0
|
66 for m = 1:c,
|
|
wolffd@0
|
67 dx = M(m,:) - x;
|
|
wolffd@0
|
68 pxm(m) = cCoeff{m} * exp(-.5 * dx * cCoinv{m} * dx');
|
|
wolffd@0
|
69 %pxm(m) = normal(dx, zeros(1,dim), diag(K(m,:)));
|
|
wolffd@0
|
70 end
|
|
wolffd@0
|
71 pxm(isnan(pxm(:))) = 0;
|
|
wolffd@0
|
72
|
|
wolffd@0
|
73 % p(x|m)
|
|
wolffd@0
|
74 if nargin>=2, Pdm(:,i) = pxm; end
|
|
wolffd@0
|
75
|
|
wolffd@0
|
76 % P(x) = P(x|M) = sum( P(m) * p(x|m) )
|
|
wolffd@0
|
77 pd(i) = P*pxm;
|
|
wolffd@0
|
78
|
|
wolffd@0
|
79 % p(m|x) = p(x|m) * P(m) / P(x)
|
|
wolffd@0
|
80 if nargout==3, pmd(:,i) = (P' .* pxm) / pd(i); end
|
|
wolffd@0
|
81
|
|
wolffd@0
|
82 end
|
|
wolffd@0
|
83
|
|
wolffd@0
|
84
|
|
wolffd@0
|
85 return;
|
|
wolffd@0
|
86
|
|
wolffd@0
|
87 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
|
wolffd@0
|
88 %
|
|
wolffd@0
|
89 % subfunction normal
|
|
wolffd@0
|
90 %
|
|
wolffd@0
|
91 % computes probability of x when mean and covariance matrix
|
|
wolffd@0
|
92 % of a distribution are known
|
|
wolffd@0
|
93
|
|
wolffd@0
|
94 function result = normal(x, mu, co)
|
|
wolffd@0
|
95
|
|
wolffd@0
|
96 [l dim] = size(x);
|
|
wolffd@0
|
97 coinv = inv(co);
|
|
wolffd@0
|
98 coeff = 1 / ((2*pi)^(dim/2)*det(co)^.5);
|
|
wolffd@0
|
99 diff = x - mu;
|
|
wolffd@0
|
100 result = coeff * exp(-.5 * diff * coinv * diff');
|
