annotate toolboxes/MIRtoolbox1.3.2/somtoolbox/som_umat.m @ 0:e9a9cd732c1e tip

first hg version after svn
author wolffd
date Tue, 10 Feb 2015 15:05:51 +0000
parents
children
rev   line source
wolffd@0 1 function U = som_umat(sMap, varargin)
wolffd@0 2
wolffd@0 3 %SOM_UMAT Compute unified distance matrix of self-organizing map.
wolffd@0 4 %
wolffd@0 5 % U = som_umat(sMap, [argID, value, ...])
wolffd@0 6 %
wolffd@0 7 % U = som_umat(sMap);
wolffd@0 8 % U = som_umat(M,sTopol,'median','mask',[1 1 0 1]);
wolffd@0 9 %
wolffd@0 10 % Input and output arguments ([]'s are optional):
wolffd@0 11 % sMap (struct) map struct or
wolffd@0 12 % (matrix) the codebook matrix of the map
wolffd@0 13 % [argID, (string) See below. The values which are unambiguous can
wolffd@0 14 % value] (varies) be given without the preceeding argID.
wolffd@0 15 %
wolffd@0 16 % U (matrix) u-matrix of the self-organizing map
wolffd@0 17 %
wolffd@0 18 % Here are the valid argument IDs and corresponding values. The values which
wolffd@0 19 % are unambiguous (marked with '*') can be given without the preceeding argID.
wolffd@0 20 % 'mask' (vector) size dim x 1, weighting factors for different
wolffd@0 21 % components (same as BMU search mask)
wolffd@0 22 % 'msize' (vector) map grid size
wolffd@0 23 % 'topol' *(struct) topology struct
wolffd@0 24 % 'som_topol','sTopol' = 'topol'
wolffd@0 25 % 'lattice' *(string) map lattice, 'hexa' or 'rect'
wolffd@0 26 % 'mode' *(string) 'min','mean','median','max', default is 'median'
wolffd@0 27 %
wolffd@0 28 % NOTE! the U-matrix is always calculated for 'sheet'-shaped map and
wolffd@0 29 % the map grid must be at most 2-dimensional.
wolffd@0 30 %
wolffd@0 31 % For more help, try 'type som_umat' or check out online documentation.
wolffd@0 32 % See also SOM_SHOW, SOM_CPLANE.
wolffd@0 33
wolffd@0 34 %%%%%%%%%%%%% DETAILED DESCRIPTION %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0 35 %
wolffd@0 36 % som_umat
wolffd@0 37 %
wolffd@0 38 % PURPOSE
wolffd@0 39 %
wolffd@0 40 % Computes the unified distance matrix of a SOM.
wolffd@0 41 %
wolffd@0 42 % SYNTAX
wolffd@0 43 %
wolffd@0 44 % U = som_umat(sM)
wolffd@0 45 % U = som_umat(...,'argID',value,...)
wolffd@0 46 % U = som_umat(...,value,...)
wolffd@0 47 %
wolffd@0 48 % DESCRIPTION
wolffd@0 49 %
wolffd@0 50 % Compute and return the unified distance matrix of a SOM.
wolffd@0 51 % For example a case of 5x1 -sized map:
wolffd@0 52 % m(1) m(2) m(3) m(4) m(5)
wolffd@0 53 % where m(i) denotes one map unit. The u-matrix is a 9x1 vector:
wolffd@0 54 % u(1) u(1,2) u(2) u(2,3) u(3) u(3,4) u(4) u(4,5) u(5)
wolffd@0 55 % where u(i,j) is the distance between map units m(i) and m(j)
wolffd@0 56 % and u(k) is the mean (or minimum, maximum or median) of the
wolffd@0 57 % surrounding values, e.g. u(3) = (u(2,3) + u(3,4))/2.
wolffd@0 58 %
wolffd@0 59 % Note that the u-matrix is always calculated for 'sheet'-shaped map and
wolffd@0 60 % the map grid must be at most 2-dimensional.
wolffd@0 61 %
wolffd@0 62 % REFERENCES
wolffd@0 63 %
wolffd@0 64 % Ultsch, A., Siemon, H.P., "Kohonen's Self-Organizing Feature Maps
wolffd@0 65 % for Exploratory Data Analysis", in Proc. of INNC'90,
wolffd@0 66 % International Neural Network Conference, Dordrecht,
wolffd@0 67 % Netherlands, 1990, pp. 305-308.
wolffd@0 68 % Kohonen, T., "Self-Organizing Map", 2nd ed., Springer-Verlag,
wolffd@0 69 % Berlin, 1995, pp. 117-119.
wolffd@0 70 % Iivarinen, J., Kohonen, T., Kangas, J., Kaski, S., "Visualizing
wolffd@0 71 % the Clusters on the Self-Organizing Map", in proceedings of
wolffd@0 72 % Conference on Artificial Intelligence Research in Finland,
wolffd@0 73 % Helsinki, Finland, 1994, pp. 122-126.
wolffd@0 74 % Kraaijveld, M.A., Mao, J., Jain, A.K., "A Nonlinear Projection
wolffd@0 75 % Method Based on Kohonen's Topology Preserving Maps", IEEE
wolffd@0 76 % Transactions on Neural Networks, vol. 6, no. 3, 1995, pp. 548-559.
wolffd@0 77 %
wolffd@0 78 % REQUIRED INPUT ARGUMENTS
wolffd@0 79 %
wolffd@0 80 % sM (struct) SOM Toolbox struct or the codebook matrix of the map.
wolffd@0 81 % (matrix) The matrix may be 3-dimensional in which case the first
wolffd@0 82 % two dimensions are taken for the map grid dimensions (msize).
wolffd@0 83 %
wolffd@0 84 % OPTIONAL INPUT ARGUMENTS
wolffd@0 85 %
wolffd@0 86 % argID (string) Argument identifier string (see below).
wolffd@0 87 % value (varies) Value for the argument (see below).
wolffd@0 88 %
wolffd@0 89 % The optional arguments are given as 'argID',value -pairs. If the
wolffd@0 90 % value is unambiguous, it can be given without the preceeding argID.
wolffd@0 91 % If an argument is given value multiple times, the last one is used.
wolffd@0 92 %
wolffd@0 93 % Below is the list of valid arguments:
wolffd@0 94 % 'mask' (vector) mask to be used in calculating
wolffd@0 95 % the interunit distances, size [dim 1]. Default is
wolffd@0 96 % the one in sM (field sM.mask) or a vector of
wolffd@0 97 % ones if only a codebook matrix was given.
wolffd@0 98 % 'topol' (struct) topology of the map. Default is the one
wolffd@0 99 % in sM (field sM.topol).
wolffd@0 100 % 'sTopol','som_topol' (struct) = 'topol'
wolffd@0 101 % 'msize' (vector) map grid dimensions
wolffd@0 102 % 'lattice' (string) map lattice 'rect' or 'hexa'
wolffd@0 103 % 'mode' (string) 'min', 'mean', 'median' or 'max'
wolffd@0 104 % Map unit value computation method. In fact,
wolffd@0 105 % eval-function is used to evaluate this, so
wolffd@0 106 % you can give other computation methods as well.
wolffd@0 107 % Default is 'median'.
wolffd@0 108 %
wolffd@0 109 % OUTPUT ARGUMENTS
wolffd@0 110 %
wolffd@0 111 % U (matrix) the unified distance matrix of the SOM
wolffd@0 112 % size 2*n1-1 x 2*n2-1, where n1 = msize(1) and n2 = msize(2)
wolffd@0 113 %
wolffd@0 114 % EXAMPLES
wolffd@0 115 %
wolffd@0 116 % U = som_umat(sM);
wolffd@0 117 % U = som_umat(sM.codebook,sM.topol,'median','mask',[1 1 0 1]);
wolffd@0 118 % U = som_umat(rand(10,10,4),'hexa','rect');
wolffd@0 119 %
wolffd@0 120 % SEE ALSO
wolffd@0 121 %
wolffd@0 122 % som_show show the selected component planes and the u-matrix
wolffd@0 123 % som_cplane draw a 2D unified distance matrix
wolffd@0 124
wolffd@0 125 % Copyright (c) 1997-2000 by the SOM toolbox programming team.
wolffd@0 126 % http://www.cis.hut.fi/projects/somtoolbox/
wolffd@0 127
wolffd@0 128 % Version 1.0beta juuso 260997
wolffd@0 129 % Version 2.0beta juuso 151199, 151299, 200900
wolffd@0 130
wolffd@0 131 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0 132 %% check arguments
wolffd@0 133
wolffd@0 134 error(nargchk(1, Inf, nargin)); % check no. of input arguments is correct
wolffd@0 135
wolffd@0 136 % sMap
wolffd@0 137 if isstruct(sMap),
wolffd@0 138 M = sMap.codebook;
wolffd@0 139 sTopol = sMap.topol;
wolffd@0 140 mask = sMap.mask;
wolffd@0 141 elseif isnumeric(sMap),
wolffd@0 142 M = sMap;
wolffd@0 143 si = size(M);
wolffd@0 144 dim = si(end);
wolffd@0 145 if length(si)>2, msize = si(1:end-1);
wolffd@0 146 else msize = [si(1) 1];
wolffd@0 147 end
wolffd@0 148 munits = prod(msize);
wolffd@0 149 sTopol = som_set('som_topol','msize',msize,'lattice','rect','shape','sheet');
wolffd@0 150 mask = ones(dim,1);
wolffd@0 151 M = reshape(M,[munits,dim]);
wolffd@0 152 end
wolffd@0 153 mode = 'median';
wolffd@0 154
wolffd@0 155 % varargin
wolffd@0 156 i=1;
wolffd@0 157 while i<=length(varargin),
wolffd@0 158 argok = 1;
wolffd@0 159 if ischar(varargin{i}),
wolffd@0 160 switch varargin{i},
wolffd@0 161 % argument IDs
wolffd@0 162 case 'mask', i=i+1; mask = varargin{i};
wolffd@0 163 case 'msize', i=i+1; sTopol.msize = varargin{i};
wolffd@0 164 case 'lattice', i=i+1; sTopol.lattice = varargin{i};
wolffd@0 165 case {'topol','som_topol','sTopol'}, i=i+1; sTopol = varargin{i};
wolffd@0 166 case 'mode', i=i+1; mode = varargin{i};
wolffd@0 167 % unambiguous values
wolffd@0 168 case {'hexa','rect'}, sTopol.lattice = varargin{i};
wolffd@0 169 case {'min','mean','median','max'}, mode = varargin{i};
wolffd@0 170 otherwise argok=0;
wolffd@0 171 end
wolffd@0 172 elseif isstruct(varargin{i}) & isfield(varargin{i},'type'),
wolffd@0 173 switch varargin{i}(1).type,
wolffd@0 174 case 'som_topol', sTopol = varargin{i};
wolffd@0 175 case 'som_map', sTopol = varargin{i}.topol;
wolffd@0 176 otherwise argok=0;
wolffd@0 177 end
wolffd@0 178 else
wolffd@0 179 argok = 0;
wolffd@0 180 end
wolffd@0 181 if ~argok,
wolffd@0 182 disp(['(som_umat) Ignoring invalid argument #' num2str(i+1)]);
wolffd@0 183 end
wolffd@0 184 i = i+1;
wolffd@0 185 end
wolffd@0 186
wolffd@0 187 % check
wolffd@0 188 [munits dim] = size(M);
wolffd@0 189 if prod(sTopol.msize)~=munits,
wolffd@0 190 error('Map grid size does not match the number of map units.')
wolffd@0 191 end
wolffd@0 192 if length(sTopol.msize)>2,
wolffd@0 193 error('Can only handle 1- and 2-dimensional map grids.')
wolffd@0 194 end
wolffd@0 195 if prod(sTopol.msize)==1,
wolffd@0 196 warning('Only one codebook vector.'); U = []; return;
wolffd@0 197 end
wolffd@0 198 if ~strcmp(sTopol.shape,'sheet'),
wolffd@0 199 disp(['The ' sTopol.shape ' shape of the map ignored. Using sheet instead.']);
wolffd@0 200 end
wolffd@0 201
wolffd@0 202 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0 203 %% initialize variables
wolffd@0 204
wolffd@0 205 y = sTopol.msize(1);
wolffd@0 206 x = sTopol.msize(2);
wolffd@0 207 lattice = sTopol.lattice;
wolffd@0 208 shape = sTopol.shape;
wolffd@0 209 M = reshape(M,[y x dim]);
wolffd@0 210
wolffd@0 211 ux = 2 * x - 1;
wolffd@0 212 uy = 2 * y - 1;
wolffd@0 213 U = zeros(uy, ux);
wolffd@0 214
wolffd@0 215 calc = sprintf('%s(a)',mode);
wolffd@0 216
wolffd@0 217 if size(mask,2)>1, mask = mask'; end
wolffd@0 218
wolffd@0 219 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0 220 %% u-matrix computation
wolffd@0 221
wolffd@0 222 % distances between map units
wolffd@0 223
wolffd@0 224 if strcmp(lattice, 'rect'), % rectangular lattice
wolffd@0 225
wolffd@0 226 for j=1:y, for i=1:x,
wolffd@0 227 if i<x,
wolffd@0 228 dx = (M(j,i,:) - M(j,i+1,:)).^2; % horizontal
wolffd@0 229 U(2*j-1,2*i) = sqrt(mask'*dx(:));
wolffd@0 230 end
wolffd@0 231 if j<y,
wolffd@0 232 dy = (M(j,i,:) - M(j+1,i,:)).^2; % vertical
wolffd@0 233 U(2*j,2*i-1) = sqrt(mask'*dy(:));
wolffd@0 234 end
wolffd@0 235 if j<y & i<x,
wolffd@0 236 dz1 = (M(j,i,:) - M(j+1,i+1,:)).^2; % diagonals
wolffd@0 237 dz2 = (M(j+1,i,:) - M(j,i+1,:)).^2;
wolffd@0 238 U(2*j,2*i) = (sqrt(mask'*dz1(:))+sqrt(mask'*dz2(:)))/(2 * sqrt(2));
wolffd@0 239 end
wolffd@0 240 end
wolffd@0 241 end
wolffd@0 242
wolffd@0 243 elseif strcmp(lattice, 'hexa') % hexagonal lattice
wolffd@0 244
wolffd@0 245 for j=1:y,
wolffd@0 246 for i=1:x,
wolffd@0 247 if i<x,
wolffd@0 248 dx = (M(j,i,:) - M(j,i+1,:)).^2; % horizontal
wolffd@0 249 U(2*j-1,2*i) = sqrt(mask'*dx(:));
wolffd@0 250 end
wolffd@0 251
wolffd@0 252 if j<y, % diagonals
wolffd@0 253 dy = (M(j,i,:) - M(j+1,i,:)).^2;
wolffd@0 254 U(2*j,2*i-1) = sqrt(mask'*dy(:));
wolffd@0 255
wolffd@0 256 if rem(j,2)==0 & i<x,
wolffd@0 257 dz= (M(j,i,:) - M(j+1,i+1,:)).^2;
wolffd@0 258 U(2*j,2*i) = sqrt(mask'*dz(:));
wolffd@0 259 elseif rem(j,2)==1 & i>1,
wolffd@0 260 dz = (M(j,i,:) - M(j+1,i-1,:)).^2;
wolffd@0 261 U(2*j,2*i-2) = sqrt(mask'*dz(:));
wolffd@0 262 end
wolffd@0 263 end
wolffd@0 264 end
wolffd@0 265 end
wolffd@0 266
wolffd@0 267 end
wolffd@0 268
wolffd@0 269 % values on the units
wolffd@0 270
wolffd@0 271 if (uy == 1 | ux == 1),
wolffd@0 272 % in 1-D case, mean is equal to median
wolffd@0 273
wolffd@0 274 ma = max([ux uy]);
wolffd@0 275 for i = 1:2:ma,
wolffd@0 276 if i>1 & i<ma,
wolffd@0 277 a = [U(i-1) U(i+1)];
wolffd@0 278 U(i) = eval(calc);
wolffd@0 279 elseif i==1, U(i) = U(i+1);
wolffd@0 280 else U(i) = U(i-1); % i==ma
wolffd@0 281 end
wolffd@0 282 end
wolffd@0 283
wolffd@0 284 elseif strcmp(lattice, 'rect')
wolffd@0 285
wolffd@0 286 for j=1:2:uy,
wolffd@0 287 for i=1:2:ux,
wolffd@0 288 if i>1 & j>1 & i<ux & j<uy, % middle part of the map
wolffd@0 289 a = [U(j,i-1) U(j,i+1) U(j-1,i) U(j+1,i)];
wolffd@0 290 elseif j==1 & i>1 & i<ux, % upper edge
wolffd@0 291 a = [U(j,i-1) U(j,i+1) U(j+1,i)];
wolffd@0 292 elseif j==uy & i>1 & i<ux, % lower edge
wolffd@0 293 a = [U(j,i-1) U(j,i+1) U(j-1,i)];
wolffd@0 294 elseif i==1 & j>1 & j<uy, % left edge
wolffd@0 295 a = [U(j,i+1) U(j-1,i) U(j+1,i)];
wolffd@0 296 elseif i==ux & j>1 & j<uy, % right edge
wolffd@0 297 a = [U(j,i-1) U(j-1,i) U(j+1,i)];
wolffd@0 298 elseif i==1 & j==1, % top left corner
wolffd@0 299 a = [U(j,i+1) U(j+1,i)];
wolffd@0 300 elseif i==ux & j==1, % top right corner
wolffd@0 301 a = [U(j,i-1) U(j+1,i)];
wolffd@0 302 elseif i==1 & j==uy, % bottom left corner
wolffd@0 303 a = [U(j,i+1) U(j-1,i)];
wolffd@0 304 elseif i==ux & j==uy, % bottom right corner
wolffd@0 305 a = [U(j,i-1) U(j-1,i)];
wolffd@0 306 else
wolffd@0 307 a = 0;
wolffd@0 308 end
wolffd@0 309 U(j,i) = eval(calc);
wolffd@0 310 end
wolffd@0 311 end
wolffd@0 312
wolffd@0 313 elseif strcmp(lattice, 'hexa')
wolffd@0 314
wolffd@0 315 for j=1:2:uy,
wolffd@0 316 for i=1:2:ux,
wolffd@0 317 if i>1 & j>1 & i<ux & j<uy, % middle part of the map
wolffd@0 318 a = [U(j,i-1) U(j,i+1)];
wolffd@0 319 if rem(j-1,4)==0, a = [a, U(j-1,i-1) U(j-1,i) U(j+1,i-1) U(j+1,i)];
wolffd@0 320 else a = [a, U(j-1,i) U(j-1,i+1) U(j+1,i) U(j+1,i+1)]; end
wolffd@0 321 elseif j==1 & i>1 & i<ux, % upper edge
wolffd@0 322 a = [U(j,i-1) U(j,i+1) U(j+1,i-1) U(j+1,i)];
wolffd@0 323 elseif j==uy & i>1 & i<ux, % lower edge
wolffd@0 324 a = [U(j,i-1) U(j,i+1)];
wolffd@0 325 if rem(j-1,4)==0, a = [a, U(j-1,i-1) U(j-1,i)];
wolffd@0 326 else a = [a, U(j-1,i) U(j-1,i+1)]; end
wolffd@0 327 elseif i==1 & j>1 & j<uy, % left edge
wolffd@0 328 a = U(j,i+1);
wolffd@0 329 if rem(j-1,4)==0, a = [a, U(j-1,i) U(j+1,i)];
wolffd@0 330 else a = [a, U(j-1,i) U(j-1,i+1) U(j+1,i) U(j+1,i+1)]; end
wolffd@0 331 elseif i==ux & j>1 & j<uy, % right edge
wolffd@0 332 a = U(j,i-1);
wolffd@0 333 if rem(j-1,4)==0, a=[a, U(j-1,i) U(j-1,i-1) U(j+1,i) U(j+1,i-1)];
wolffd@0 334 else a = [a, U(j-1,i) U(j+1,i)]; end
wolffd@0 335 elseif i==1 & j==1, % top left corner
wolffd@0 336 a = [U(j,i+1) U(j+1,i)];
wolffd@0 337 elseif i==ux & j==1, % top right corner
wolffd@0 338 a = [U(j,i-1) U(j+1,i-1) U(j+1,i)];
wolffd@0 339 elseif i==1 & j==uy, % bottom left corner
wolffd@0 340 if rem(j-1,4)==0, a = [U(j,i+1) U(j-1,i)];
wolffd@0 341 else a = [U(j,i+1) U(j-1,i) U(j-1,i+1)]; end
wolffd@0 342 elseif i==ux & j==uy, % bottom right corner
wolffd@0 343 if rem(j-1,4)==0, a = [U(j,i-1) U(j-1,i) U(j-1,i-1)];
wolffd@0 344 else a = [U(j,i-1) U(j-1,i)]; end
wolffd@0 345 else
wolffd@0 346 a=0;
wolffd@0 347 end
wolffd@0 348 U(j,i) = eval(calc);
wolffd@0 349 end
wolffd@0 350 end
wolffd@0 351 end
wolffd@0 352
wolffd@0 353 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0 354 %% normalization between [0,1]
wolffd@0 355
wolffd@0 356 % U = U - min(min(U));
wolffd@0 357 % ma = max(max(U)); if ma > 0, U = U / ma; end
wolffd@0 358
wolffd@0 359 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0 360
wolffd@0 361
wolffd@0 362