wolffd@0
|
1 function [K,P] = som_estimate_gmm(sM, sD)
|
wolffd@0
|
2
|
wolffd@0
|
3 %SOM_ESTIMATE_GMM Estimate a gaussian mixture model based on map.
|
wolffd@0
|
4 %
|
wolffd@0
|
5 % [K,P] = som_estimate_gmm(sM, sD)
|
wolffd@0
|
6 %
|
wolffd@0
|
7 % Input and output arguments:
|
wolffd@0
|
8 % sM (struct) map struct
|
wolffd@0
|
9 % sD (struct) data struct
|
wolffd@0
|
10 % (matrix) size dlen x dim, the data to use when estimating
|
wolffd@0
|
11 % the gaussian kernels
|
wolffd@0
|
12 %
|
wolffd@0
|
13 % K (matrix) size munits x dim, kernel width parametes for
|
wolffd@0
|
14 % each map unit
|
wolffd@0
|
15 % P (vector) size 1 x munits, a priori probability of each map unit
|
wolffd@0
|
16 %
|
wolffd@0
|
17 % See also SOM_PROBABILITY_GMM.
|
wolffd@0
|
18
|
wolffd@0
|
19 % Reference: Alhoniemi, E., Himberg, J., Vesanto, J.,
|
wolffd@0
|
20 % "Probabilistic measures for responses of Self-Organizing Maps",
|
wolffd@0
|
21 % Proceedings of Computational Intelligence Methods and
|
wolffd@0
|
22 % Applications (CIMA), 1999, Rochester, N.Y., USA, pp. 286-289.
|
wolffd@0
|
23
|
wolffd@0
|
24 % Contributed to SOM Toolbox vs2, February 2nd, 2000 by Esa Alhoniemi
|
wolffd@0
|
25 % Copyright (c) by Esa Alhoniemi
|
wolffd@0
|
26 % http://www.cis.hut.fi/projects/somtoolbox/
|
wolffd@0
|
27
|
wolffd@0
|
28 % ecco 180298 juuso 050100 250400
|
wolffd@0
|
29
|
wolffd@0
|
30 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
wolffd@0
|
31
|
wolffd@0
|
32 [c, dim] = size(sM.codebook);
|
wolffd@0
|
33 M = sM.codebook;
|
wolffd@0
|
34
|
wolffd@0
|
35 if isstruct(sD), D = sD.data; else D = sD; end
|
wolffd@0
|
36 dlen = length(D(:,1));
|
wolffd@0
|
37
|
wolffd@0
|
38 %%%%%%%%%%%%%%%%%%%%%
|
wolffd@0
|
39 % compute hits & bmus
|
wolffd@0
|
40
|
wolffd@0
|
41 [bmus, qerrs] = som_bmus(sM, D);
|
wolffd@0
|
42 hits = zeros(1,c);
|
wolffd@0
|
43 for i = 1:c, hits(i) = sum(bmus == i); end
|
wolffd@0
|
44
|
wolffd@0
|
45 %%%%%%%%%%%%%%%%%%%%
|
wolffd@0
|
46 % a priori
|
wolffd@0
|
47
|
wolffd@0
|
48 % neighborhood kernel
|
wolffd@0
|
49 r = sM.trainhist(end).radius_fin; % neighborhood radius
|
wolffd@0
|
50 if isempty(r) | isnan(r), r=1; end
|
wolffd@0
|
51 Ud = som_unit_dists(sM);
|
wolffd@0
|
52 Ud = Ud.^2;
|
wolffd@0
|
53 r = r^2;
|
wolffd@0
|
54 if r==0, r=eps; end % to get rid of div-by-zero errors
|
wolffd@0
|
55 switch sM.neigh,
|
wolffd@0
|
56 case 'bubble', H = (Ud<=r);
|
wolffd@0
|
57 case 'gaussian', H = exp(-Ud/(2*r));
|
wolffd@0
|
58 case 'cutgauss', H = exp(-Ud/(2*r)) .* (Ud<=r);
|
wolffd@0
|
59 case 'ep', H = (1-Ud/r) .* (Ud<=r);
|
wolffd@0
|
60 end
|
wolffd@0
|
61
|
wolffd@0
|
62 % a priori prob. = hit histogram weighted by the neighborhood kernel
|
wolffd@0
|
63 P = hits*H;
|
wolffd@0
|
64 P = P/sum(P);
|
wolffd@0
|
65
|
wolffd@0
|
66 %%%%%%%%%%%%%%%%%%%%
|
wolffd@0
|
67 % kernel widths (& centers)
|
wolffd@0
|
68
|
wolffd@0
|
69 K = ones(c, dim) * NaN; % kernel widths
|
wolffd@0
|
70 for m = 1:c,
|
wolffd@0
|
71 w = H(bmus,m);
|
wolffd@0
|
72 w = w/sum(w);
|
wolffd@0
|
73 for i = 1:dim,
|
wolffd@0
|
74 d = (D(:,i) - M(m,i)).^2; % compute variance of ith
|
wolffd@0
|
75 K(m,i) = w'*d; % variable of centroid m
|
wolffd@0
|
76 end
|
wolffd@0
|
77 end
|
wolffd@0
|
78
|
wolffd@0
|
79 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|