Reading material¶
[KGB06]: first fMRI searchlight paper. [OWD+10]: surface based searchlight (outside the scope of this exercise).
Write a searchlight function that computes a generic dataset measure¶
Volume based searchlight analysis proceeds by defining a local neighborhood of voxels around each voxel in the brain volume. This subset of voxels can be thought of as a mini dataset, where the mask that defines the voxels included in the dataset is the searchlight sphere. Because we can treat each searchlight as a dataset, we can build a searchlight function that will compute any dataset measure that we specify. This allows us to reuse code, and run searchlights for different purposes.
We have provided a couple of helper functions that do some of the heavy lifting involved defining the sets of neighborhood voxels. For this, use cosmo spherical neighborhood that was presented in the previous exercise.
With the help of these functions, write a generic searchlight function that satisfies the following definition:
function results_map = cosmo_searchlight(ds, nbrhood, measure, varargin)
% Generic searchlight function returns a map of results computed at each
% searchlight location
%
% results_map=cosmo_searchlight(dataset, nbrhood, measure, ...)
%
% Inputs:
% ds dataset struct with field .samples (NSxNF)
% nbrhood Neighborhood structure with fields:
% .a struct with dataset attributes
% .fa struct with feature attributes. Each field
% should have NF values in the second dimension
% .neighbors cell with NF mappings from center_ids in output
% dataset to feature ids in input dataset.
% Suitable neighborhood structs can be generated
% using:
% - cosmo_spherical_neighborhood (fmri volume)
% - cosmo_surficial_neighborhood (fmri surface)
% - cosmo_meeg_chan_neigborhood (MEEG channels)
% - cosmo_interval_neighborhood (MEEG time, freq)
% - cosmo_cross_neighborhood (to cross neighborhoods
% from the neighborhood
% functions above)
% measure function handle to a dataset measure. A dataset
% measure has the function signature:
% output = measure(dataset, args)
% where output must be a struct with fields .samples
% (as a column vector) and optionally a field .sa
% with sample attributes.
% Typical measures are:
% - cosmo_correlation_measure
% - cosmo_crossvalidation_measure
% - cosmo_target_dsm_corr_measure
% 'center_ids', ids vector indicating center ids to be used as a
% searchlight center. By default all feature ids are
% used (i.e. ids=1:numel(nbrhood.neighbors). The
% output results_map.samples has size N in the 2nd
% dimension.
% 'progress', p Show progress every p steps
% 'nproc', np If the Matlab parallel processing toolbox, or the
% GNU Octave parallel package is available, use
% np parallel threads. (Multiple threads may speed
% up searchlight computations).
% If parallel processing is not available, or if
% this option is not provided, then a single thread
% is used.
% K, V any key-value pair (K,V) with arguments for the
% measure function handle. Alternatively a struct
% can be used
%
% Output:
% results_map a dataset struct where the samples
% contain the results of the searchlight analysis.
% If measure returns datasets all of size Nx1 and
% there are M center_ids
% (M=numel(nbrhood.neighbors) if center_ids is not
% provided), then results_map.samples has size MxN.
% If nbrhood has fields .a and .fa, these are part
% of the output (with .fa sliced according to
% center_ids)
%
% Example:
% % use a minimal dataset with 6 voxels
% ds=cosmo_synthetic_dataset('nchunks',5);
% %
% % define neighborhood (progress is set to false to suppress output)
% radius=1; % radius=3 is typical for fMRI datasets
% nbrhood=cosmo_spherical_neighborhood(ds,'radius',radius,...
% 'progress',false);
% %
% % define measure and its arguments; here crossvalidation with LDA
% % classifier to compute classification accuracies
% args=struct();
% args.classifier = @cosmo_classify_lda;
% args.partitions = cosmo_nfold_partitioner(ds);
% measure=@cosmo_crossvalidation_measure;
% %
% % run searchlight (without showing progress bar)
% result=cosmo_searchlight(ds,nbrhood,measure,'progress',0,args);
% %
% % show results:
% % - .samples contains classification accuracy
% % - .fa.nvoxels is the number of voxels in each searchlight
% % - .fa.radius is the radius of each searchlight
% cosmo_disp(result)
% %|| .a
% %|| .fdim
% %|| .labels
% %|| { 'i' 'j' 'k' }
% %|| .values
% %|| { [ 1 2 3 ] [ 1 2 ] [ 1 ] }
% %|| .vol
% %|| .mat
% %|| [ 2 0 0 -3
% %|| 0 2 0 -3
% %|| 0 0 2 -3
% %|| 0 0 0 1 ]
% %|| .dim
% %|| [ 3 2 1 ]
% %|| .xform
% %|| 'scanner_anat'
% %|| .fa
% %|| .nvoxels
% %|| [ 3 4 3 3 4 3 ]
% %|| .radius
% %|| [ 1 1 1 1 1 1 ]
% %|| .center_ids
% %|| [ 1 2 3 4 5 6 ]
% %|| .i
% %|| [ 1 2 3 1 2 3 ]
% %|| .j
% %|| [ 1 1 1 2 2 2 ]
% %|| .k
% %|| [ 1 1 1 1 1 1 ]
% %|| .samples
% %|| [ 1 1 1 0.9 1 0.7 ]
% %|| .sa
% %|| .labels
% %|| { 'accuracy' }
%
% Notes:
% - neighborhoods can be defined using one or more of the
% cosmo_*_neighborhood functions
%
% See also: cosmo_correlation_measure,
% cosmo_crossvalidation_measure,
% cosmo_dissimilarity_matrix_measure,
% cosmo_spherical_neighborhood,cosmo_surficial_neighborhood,
% cosmo_meeg_chan_neigborhood, cosmo_interval_neighborhood
% cosmo_cross_neighborhood
%
% # For CoSMoMVPA's copyright information and license terms, #
% # see the COPYING file distributed with CoSMoMVPA. #
Hint: cosmo searchlight skl
Solution: cosmo searchlight
