cosmo find local extrema hdrΒΆ

function [feature_ids,scores]=cosmo_find_local_extrema(ds, nbrhood, varargin)
% find local extrema in a dataset using a neighborhood
%
% [feature_ids,scores]=cosmo_find_local_extrema(ds, nbrhood, ...)
%
% Inputs:
%   ds              dataset struct with one sample, i.e. size(ds.samples,1)
%                   must be 1. NaN values ds.samples are ignored.
%   nbrhood         neighborhood structure corresponding with dataset, with
%                   a field .neighbors so that .neighbors{k} contains the
%                   feature ids of neighbors of feature k.
%   'fitness', f    optional function handle so that [score,id]=f(x)
%                   returns the index with the maximum 'feature score'.
%                   The default is @(x)max(x,[],2), i.e. a function that
%                   returns the maximum value and the index of the maximum
%                   value of the input.
%   'count', c      optional number of elements to select. The default is
%                   Inf.
%
% Output:
%   feature_ids     1xN vector with feature ids.
%   scores          1xN vector with fitness scores. It holds for all I and
%                   J in 1:N that, if I<J, then:
%                   - fitness(ds.samples(I))>fitness(ds.samples(J))
%                   - intersect(nbrhood.neighbors{I},nbrhood.neighbors{J})
%                       is empty
%                   - if I_N=setdiff(nbrhood.neighbors{I},I), then
%                     fitness(ds.samples(I_N))<=fitness(ds.samples(I))
%
%
% Example:
%     % Generate tiny dataset with 6 voxels
%     ds=cosmo_synthetic_dataset('ntargets',1,'nchunks',1);
%     cosmo_disp(ds.samples)
%     %|| [ 2.03     0.584     -1.44    -0.518      1.19     -1.33 ]
%     %
%     % show unflattened shape
%     cosmo_disp(squeeze(cosmo_unflatten(ds)))
%     %|| [  2.03    -0.518
%     %||   0.584      1.19
%     %||   -1.44     -1.33 ]
%     %
%     nh=cosmo_spherical_neighborhood(ds,'radius',1,'progress',false);
%     % find local maxima within neighborhood of 1 voxel radius
%     [feature_ids,scores]=cosmo_find_local_extrema(ds,nh);
%     cosmo_disp(feature_ids)
%     %|| [ 1         5         3 ]
%     cosmo_disp(scores)
%     %|| [ 2.03      1.19     -1.44 ]
%     %
%     % only return two feature ids
%     [feature_ids,scores]=cosmo_find_local_extrema(ds,nh,'count',2);
%     cosmo_disp(feature_ids)
%     %|| [ 1         5 ]
%     cosmo_disp(scores)
%     %|| [ 2.03      1.19 ]
%     %
%     % use another fitness function, namely local minima
%     [feature_ids,scores]=cosmo_find_local_extrema(ds,nh,'fitness',@min);
%     cosmo_disp(feature_ids)
%     %|| [ 3         4 ]
%     cosmo_disp(scores)
%     %|| [ -1.44    -0.518 ]
%     %
%     % find local maxima within neighborhood of 2 voxel radius
%     nh=cosmo_spherical_neighborhood(ds,'radius',2,'progress',false);
%     [feature_ids,scores]=cosmo_find_local_extrema(ds,nh);
%     cosmo_disp(feature_ids)
%     %|| [ 1         6 ]
%     cosmo_disp(scores)
%     %|| [ 2.03     -1.33 ]
%
%
%     % A more sophisticated example showing how a table of local maxima
%     % can be produced:
%     %
%     % Generate synthetic dataset with 700 voxels
%     ds=cosmo_synthetic_dataset('size','big',...
%                                        'ntargets',1,'nchunks',1);
%     %
%     % ensure LPI orientation
%     ds=cosmo_fmri_reorient(ds,'LPI');
%     %
%     % set minimum distance (in voxels) between local extrama (peaks)
%     min_voxels_between_extrema=5;
%     nh=cosmo_spherical_neighborhood(ds,...
%                         'radius',min_voxels_between_extrema,...
%                         'progress',false);
%     %
%     % find local maxima
%     % (note: local minima are not found; this would require a
%     %       'fitness' argument different from the default)
%     [feature_ids,scores]=cosmo_find_local_extrema(ds,nh);
%     %
%     % define and apply threshold
%     thr=1.96;
%     msk=scores>thr;
%     feature_ids=feature_ids(msk);
%     scores=scores(msk);
%     %
%     % get coordinates of local extrama
%     xyz=cosmo_vol_coordinates(ds,feature_ids);
%     %
%     % show table with coordinates and value at local extrema
%     fprintf('x=% 3d y=% 3d z=% 3d   score=%.3f\n',[xyz; scores]);
%     %|| x= 11 y=  7 z=  3   score=3.679
%     %|| x= -1 y=  3 z= -1   score=3.458
%     %|| x= 27 y=  7 z=  7   score=3.361
%     %|| x= 35 y= -1 z=  3   score=3.288
%     %|| x= 23 y=  1 z= -1   score=2.899
%     %|| x= 35 y= 11 z= -1   score=2.268
%
% Notes:
%   - this function can be used to define regions of interest in a
%     reproducible and objective manner.
%   - to ignore particular features, set their value to NaN
%
% #   For CoSMoMVPA's copyright information and license terms,   #
% #   see the COPYING file distributed with CoSMoMVPA.           #