run compare dsm skl

%% RSA Tutorial
% Compare DSMs
%
% #   For CoSMoMVPA's copyright information and license terms,   #
% #   see the COPYING file distributed with CoSMoMVPA.           #

subjects = {'s01', 's02', 's03', 's04', 's05', 's06', 's07', 's08'};
masks = {'ev_mask.nii', 'vt_mask.nii'};

config = cosmo_config();
study_path = fullfile(config.tutorial_data_path, 'ak6');

%%
% In a nested loop over masks then subjects:
%       load each dataset
%       demean it, get DSM, save it in dsms

n_subjects = numel(subjects);
n_masks = numel(masks);

counter = 0;

for m = 1:n_masks
    msk = masks{m};
    for s = 1:length(subjects)
        sub = subjects{s};
        sub_path = fullfile(study_path, sub);

        % load dataset
        ds_fn = fullfile(sub_path, 'glm_T_stats_perrun.nii');
        mask_fn = fullfile(sub_path, msk);
        ds_full = cosmo_fmri_dataset(ds_fn, ...
                                     'mask', mask_fn, ...
                                     'targets', repmat(1:6, 1, 10)');

        % compute average for each unique target
        ds = cosmo_fx(ds_full, @(x)mean(x, 1), 'targets', 1);

        % remove constant features
        ds = cosmo_remove_useless_data(ds);

        % demean
        % Comment this out to see the effects of demeaning vs. not
        ds.samples = bsxfun(@minus, ds.samples, mean(ds.samples, 1));

        % compute the one-minus-correlation value for each pair of
        % targets.
        % (Hint: use cosmo_pdist with the 'correlation' argument)
        %%%% >>> Your code here <<< %%%%

        if counter == 0
            % first dsm, allocate space
            n_pairs = numel(dsm);
            neural_dsms = zeros(n_subjects * n_masks, n_pairs);
        end

        % increase counter and store the dsm as the counter-th row in
        % 'neural_dsms'
        %%%% >>> Your code here <<< %%%%
    end
end

%%
% Then add the v1 model and behavioral DSMs
models_path = fullfile(study_path, 'models');
load(fullfile(models_path, 'v1_model.mat'));
load(fullfile(models_path, 'behav_sim.mat'));
% add to dsms (hint: use comso_squareform)
%%%% >>> Your code here <<< %%%%

%%
% Now visualize the cross-correlation matrix. Remember that 'cosmo_corr'
% (or the builtin 'corr') calculates correlation coefficients between
% columns and we want between rows, so the data has to be transposed.

%%%% >>> Your code here <<< %%%%

%%
% Now use the values in the last two rows of the cross correlation matrix to
% visualize the distributions in correlations between the neural similarities
% and the v1 model/behavioral ratings. Store the result in a matrix
% 'cc_models', which should have 8 rows (corresponding to the participants)
% and 4 columns (corresponding to EV and VT correlated with both models)
%
% Rows  1 to  8: EV neural similarities
% Rows  9 to 16: VT neural similarities
% Row        17: EV model
% Row        18: behavioural similarities

%%%% >>> Your code here <<< %%%%
labels = {'v1 model~EV', 'v1 model~VT', 'behav~EV', 'behav~VT'};
figure();
boxplot(cc_models);
set(gca, 'XTick', [1:4], 'XTickLabel', labels);