spm_compact_svd.m 1.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051
  1. function U = spm_compact_svd(Y,xyz,nu)
  2. % local SVD with compact support for large matrices
  3. % FORMAT U = spm_compact_svd(Y,xyz,nu)
  4. % Y - matrix
  5. % xyz - location
  6. % nu - number of vectors
  7. %__________________________________________________________________________
  8. % Copyright (C) 2008 Wellcome Trust Centre for Neuroimaging
  9. % Karl Friston
  10. % $Id: spm_compact_svd.m 5219 2013-01-29 17:07:07Z spm $
  11. % get orders
  12. %--------------------------------------------------------------------------
  13. ns = size(Y,1); % number of samples
  14. nv = size(Y,2); % number of components (voxels/vertices)
  15. % get kernel (compact vectors)
  16. %--------------------------------------------------------------------------
  17. nc = max(fix(nv/nu),1); % voxels in compact support
  18. C = sum(Y.^2); % variance of Y
  19. U = spalloc(nv,nu,nc*nu);
  20. J = 1:nv;
  21. for i = 1:nu
  22. % find maximum variance voxel
  23. %----------------------------------------------------------------------
  24. [v,j] = max(C);
  25. d = 0;
  26. for k = 1:size(xyz,1)
  27. d = d + (xyz(k,:) - xyz(k,j)).^2;
  28. end
  29. [d,j] = sort(d);
  30. try
  31. j = j(1:nc);
  32. end
  33. % save principal eigenvector
  34. %----------------------------------------------------------------------
  35. k = J(j);
  36. u = spm_svd(Y(:,k)');
  37. U(k,i) = u(:,1);
  38. % remove compact support voxels and start again
  39. %----------------------------------------------------------------------
  40. J(j) = [];
  41. C(j) = [];
  42. xyz(:,j) = [];
  43. end