Browse code

starting som prediction fine-tuned class-performance visualisation

git-svn-id: https://svn.discofish.de/MATLAB/spmtoolbox/SVMCrossVal@112 83ab2cfd-5345-466c-8aeb-2b2739fb922d

Christoph Budziszewski authored on21/01/2009 16:34:25
Showing1 changed files
1 1
new file mode 100644
... ...
@@ -0,0 +1,207 @@
1
+function [Err,sPropTotal,sPropMunits,sPropComps] = som_distortion3(sM,D,rad)
2
+
3
+%SOM_DISTORTION3 Map distortion measures.
4
+%  
5
+% [sE,Err] = som_distortion3(sM,[D],[rad]);
6
+% 
7
+%  sE = som_distortion3(sM); 
8
+%
9
+%  Input and output arguments ([]'s are optional): 
10
+%   sM          (struct) map struct
11
+%   [D]         (matrix) a matrix, size dlen x dim
12
+%               (struct) data or map struct
13
+%                        by default the map struct is used
14
+%   [rad]       (scalar) neighborhood radius, looked from sM.trainhist
15
+%                        by default, or = 1 if that has no valid values
16
+%                           
17
+%   Err         (matrix) size munits x dim x 3
18
+%                        distortion error elements (quantization error, 
19
+%                        neighborhood bias, and neighborhood variance)
20
+%                        for each map unit and component
21
+%   sPropTotal  (struct) .n   = length of data
22
+%                        .h   = mean neighborhood function value
23
+%                        .err = errors
24
+%   sPropMunits (struct) .Ni  = hits per map unit
25
+%                        .Hi  = sum of neighborhood values for each map unit
26
+%                        .Err = errors per map unit
27
+%   sPropComps  (struct) .e1  = total squared distance to centroid
28
+%                        .eq  = total squared distance to BMU
29
+%                        .Err = errors per component
30
+%
31
+% See also  SOM_QUALITY.
32
+
33
+% Contributed to SOM Toolbox 2.0, January 3rd, 2002 by Juha Vesanto
34
+% Copyright (c) by Juha Vesanto
35
+% http://www.cis.hut.fi/projects/somtoolbox/
36
+
37
+% Version 2.0beta juuso 030102
38
+
39
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
40
+%% arguments
41
+
42
+% map
43
+[munits dim] = size(sM.codebook);
44
+
45
+% neighborhood radius
46
+if nargin<3, 
47
+  if ~isempty(sM.trainhist), 
48
+    rad = sM.trainhist(end).radius_fin; 
49
+  else 
50
+    rad = 1; 
51
+  end
52
+end
53
+if rad<eps, rad = eps; end
54
+if isempty(rad) | isnan(rad), rad = 1; end
55
+
56
+% neighborhood function
57
+Ud = som_unit_dists(sM.topol); 
58
+switch sM.neigh, 
59
+ case 'bubble',   H = (Ud <= rad);
60
+ case 'gaussian', H = exp(-(Ud.^2)/(2*rad*rad)); 
61
+ case 'cutgauss', H = exp(-(Ud.^2)/(2*rad*rad)) .* (Ud <= rad);
62
+ case 'ep',       H = (1 - (Ud.^2)/rad) .* (Ud <= rad);
63
+end  
64
+Hi = sum(H,2); 
65
+
66
+% data
67
+if nargin<2, D = sM.codebook; end
68
+if isstruct(D), D = D.data; end
69
+[dlen dim] = size(D);
70
+
71
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
72
+%% quality measures
73
+
74
+% find Voronoi sets, and calculate their properties
75
+
76
+[bmus,qerr] = som_bmus(sM,D); 
77
+M  = sM.codebook; 
78
+Vn = M; 
79
+Vm = M; 
80
+Ni = zeros(munits,dim);
81
+for i=1:munits, 
82
+  inds    = find(bmus==i);   
83
+  Ni(i,:) = sum(isfinite(D(inds,:)),1);                      % size of Voronoi set
84
+  if any(Ni(i,:)), Vn(i,:) = centroid(D(inds,:),M(i,:)); end % centroid of Voronoi set  
85
+  Vm(i,:) = centroid(M,M(i,:),H(i,:)');                      % centroid of neighborhood
86
+end
87
+
88
+HN = repmat(Hi,1,dim).*Ni; 
89
+
90
+%% distortion
91
+
92
+% quantization error (in each Voronoi set and for each component)
93
+
94
+Eqx           = zeros(munits,dim); 
95
+Dx            = (Vn(bmus,:) - D).^2; 
96
+Dx(isnan(Dx)) = 0; 
97
+for i = 1:dim, 
98
+  Eqx(:,i)    = full(sum(sparse(bmus,1:dlen,Dx(:,i),munits,dlen),2)); 
99
+end
100
+Eqx           = repmat(Hi,1,dim).*Eqx; 
101
+  
102
+% bias in neighborhood (in each Voronoi set / component)
103
+
104
+Enb = (Vn-Vm).^2;
105
+Enb = HN.*Enb; 
106
+
107
+% variance in neighborhood (in each Voronoi set / component)
108
+
109
+Env = zeros(munits,dim);
110
+for i=1:munits, Env(i,:) = H(i,:)*(M-Vm(i*ones(munits,1),:)).^2; end
111
+Env = Ni.*Env; 
112
+
113
+% total distortion (in each Voronoi set / component)
114
+
115
+Ed = Eqx + Enb + Env;
116
+
117
+%% other error measures
118
+
119
+% squared quantization error (to data centroid)
120
+
121
+me            = centroid(D,mean(M));
122
+Dx            = D - me(ones(dlen,1),:); 
123
+Dx(isnan(Dx)) = 0; 
124
+e1            = sum(Dx.^2,1); 
125
+
126
+% squared quantization error (to map units)
127
+
128
+Dx            = D - M(bmus,:);
129
+Dx(isnan(Dx)) = 0; 
130
+eq            = sum(Dx.^2,1);
131
+
132
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
133
+%% output
134
+
135
+% distortion error matrix
136
+
137
+Err  = zeros(munits,dim,5); 
138
+Err(:,:,1) = Eqx; 
139
+Err(:,:,2) = Enb; 
140
+Err(:,:,3) = Env; 
141
+
142
+% total errors
143
+
144
+sPropTotal = struct('n',sum(Ni),'h',mean(Hi),'e1',sum(e1),'err',sum(sum(Err,2),1));
145
+
146
+% properties of map units
147
+
148
+sPropMunits = struct('Ni',[],'Hi',[],'Err',[]); 
149
+sPropMunits.Ni  = Ni; 
150
+sPropMunits.Hi  = Hi; 
151
+sPropMunits.Err = squeeze(sum(Err,2));
152
+
153
+% properties of components
154
+
155
+sPropComps = struct('Err',[],'e1',[],'eq',[]);
156
+sPropComps.Err = squeeze(sum(Err,1));
157
+sPropComps.e1  = e1; 
158
+sPropComps.eq  = eq;
159
+
160
+
161
+return; 
162
+
163
+
164
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%55
165
+%% subfunctions
166
+
167
+function v = centroid(D,default,weights)
168
+  
169
+  [n dim] = size(D);
170
+  I       = sparse(isnan(D));
171
+  D(I)    = 0;
172
+  
173
+  if nargin==3, 
174
+    W    = weights(:,ones(1,dim)); 
175
+    W(I) = 0; 
176
+    D    = D.*W;
177
+    nn   = sum(W,1);
178
+  else
179
+    nn   = n-sum(I,1);
180
+  end 
181
+
182
+  c    = sum(D,1);
183
+  v    = default; 
184
+  i    = find(nn>0); 
185
+  v(i) = c(i)./nn(i);
186
+      
187
+  return; 
188
+
189
+
190
+function vis
191
+
192
+  figure
193
+  som_show(sM,'color',{Hi,'Hi'},'color',{Ni,'hits'},...
194
+           'color',{Ed,'distortion'},'color',{Eqx,'qxerror'},...
195
+           'color',{Enb,'N-bias'},'color',{Env,'N-Var'});
196
+
197
+  ed = Eqx + Enb + Env;
198
+  i = find(ed>0); 
199
+  eqx = 0*ed; eqx(i) = Eqx(i)./ed(i);
200
+  enb = 0*ed; enb(i) = Enb(i)./ed(i);
201
+  env = 0*ed; env(i) = Env(i)./ed(i);
202
+
203
+  figure
204
+  som_show(sM,'color',Hi,'color',Ni,'color',Ed,...
205
+           'color',eqx,'color',enb,'color',env); 
206
+
207
+