%% Compute first-level Item-Specific Reinstatement X Condition statistics
% specify :
%       subject (number)
%       width of rsa wins (e.g. 400), steps of rsa (e.g. 50) -> num input!
%       control ('category','memory','catmem','all'), selects DIFFs
%       diff_type ( "row", "col", "cross") -> DIFF
%       atlas Type ('AAL' or 'BA') and ROI ('Early Visual' (AAL) or 'Early Visual V2' (BA))
%       number of permutation repetitions
% run as e.g. ERS_first_level_stats(1,200,50,1500,'all','row','AAL','Early Visual',1000)
%
% computes tests of memory, category, complexity, memorability on
% item-specific scores!
% does permuation in place, each time shuffling item-specific vector
%
% computes empirical first-level statistics, several measures -> SAME,
% DIFF, Item-Specific

%SAME: Rem vs Forg, In vs Out, High Comp vs Low Comp, High Mem vs Low Mem
%DIFF (off diags): Memory, Category, Complexity, Memorability
%Item-Specific (col, row, cross), all_diff, memcat, imgcat, both
% again Memory, Category, Complexity Memorability

%item-specific metrics: On - off -> specify, row, col / cross
%control condition: off (same category, memory, both, none)

%empirical test: on-off remembered vs forgotten
%permutation: shuffle on-off values between memory categories

function ERS_first_level_itemspecific_memory(subject,control,diff_type,nrep,rsa_tag,roi_tag,stats_tag)
    %diff type: which off-items are selected -> "row", "col", "cross"
    %control: off-items selection -> "category","memory","catmem","all"
    %categories
    
    %% set basic stuff
    fs = '\';
    dataDir = '..\Data';
    ftpath = 'C:\Program Files\MATLAB\R2020b\toolbox\fieldtrip';
    save_to_disk = 1;
     
    %load folder names of subjects
    load('..\subject_list.mat'); %subNames for easier folder accessing
    
    %find correct folder of RSA data based on function input
    fileName = subNames{subject};
    disp(fileName)
    
    %% load ERS data
    ersDir = ['stats' fs 'ERS' fs rsa_tag fs roi_tag];
            
    %if subject contributes, continue with loading ERS data and info 
    cd([dataDir fs fileName fs ersDir])
    load ERSinfo.mat    
    %load ERS data and reduce to current ROI
    load ERS_tempGen.mat

    %% get indices of trial-specs (from info.trialinfo (from encRSA))
    remembered_ix = find(info.trialinfo(:,2) == 1);
    forgotten_ix = find(info.trialinfo(:,2) == 2);
    indoor_ix = find(info.trialinfo(:,5) == 1);
    outdoor_ix = find(info.trialinfo(:,5) == 2);
    
    %% Prepare template matrix for proper trial indexing (SAME - DIFF)
    tempMat = zeros(size(ERS,3)); %trialXtrial

    if strcmp(control,'category') % same - diff from same image category
            tempMat(indoor_ix,indoor_ix) = 1;
            tempMat(outdoor_ix,outdoor_ix) = 1; 
    elseif strcmp(control,'memory')   % same - diff from same memory category
            tempMat(remembered_ix,remembered_ix) = 1;
            tempMat(forgotten_ix,forgotten_ix) = 1;
    elseif strcmp(control,'memcat')   % same - diff from same image+memory category
            tempMat(intersect(indoor_ix,remembered_ix),intersect(indoor_ix,remembered_ix)) = 1;
            tempMat(intersect(indoor_ix,forgotten_ix),intersect(indoor_ix,forgotten_ix)) = 1;
            tempMat(intersect(outdoor_ix,remembered_ix),intersect(outdoor_ix,remembered_ix)) = 1;
            tempMat(intersect(outdoor_ix,forgotten_ix),intersect(outdoor_ix,forgotten_ix)) = 1;
    elseif strcmp(control,'all')      % same - all diff
            tempMat = ones(size(ERS,3));
    end
    
    %now create indices of item-specific "offs" depending just created
    %tempMat 
    off_idxs = {};

    for same = 1:size(ERS,3)
        offs = find(tempMat(same,:));  %includes SAME item index
        offs(offs==same) = [];               % dropped SAME item index
        off_idxs{same} = offs;
    end
    
    %% Calculate Item-Specific Reinstatement score
    %SAME - DIFF
    
    %store same value, avg of diffs and item-score (SAME-DIFF)
    Reinstatement = {};
    Reinstatement.SAME = zeros(size(ERS,1),size(ERS,2),size(ERS,3));
    Reinstatement.DIFF = zeros(size(ERS,1),size(ERS,2),size(ERS,3));
    Reinstatement.ITEM = zeros(size(ERS,1),size(ERS,2),size(ERS,3));
    
    for win1 = 1:size(ERS,1) %idx1 is now time
        for win2 = 1:size(ERS,2)
            for stim = 1:size(ERS,3)
                on_element = squeeze(ERS(win1,win2,stim,stim));

                %calc off elements as avg over row, col or crossed
                if strcmp(diff_type,'row')
                    off_elements = nanmean(ERS(win1,win2,stim,off_idxs{stim}));
                elseif strcmp(diff_type,'col')
                    off_elements = nanmean(ERS(win1,win2,off_idxs{stim},stim));
                elseif strcmp(diff_type,'cross')
                    off_elementsrow = squeeze(ERS(win1,win2,stim,off_idxs{stim}));
                    off_elementscol = squeeze(ERS(win1,win2,off_idxs{stim},stim));
                    off_elements = nanmean([off_elementsrow; off_elementscol]);
                end
                
                Reinstatement.SAME(win1,win2,stim) = on_element; 
                Reinstatement.DIFF(win1,win2,stim) = off_elements; %avg of all offs!
                Reinstatement.ITEM(win1,win2,stim) = on_element - off_elements;
            end
        end
    end
    
    %% Now perform multiple tests on the item-specific reinstatement
    %%% values
    % memory (rem vs forg), category (in vs out), complexity (high vs low),
    % memorability (high vs low)
    % also compute permutation stats in place, each time shuffling the
    % array

    %SameDiff diag holds on - avg(off) values. calc
    %ttest2(correct,incorrect)
    %emp_stats = zeros(size(ERS,1),5); % 3 columns: same, diff, h, p, t
    stats = {};
    stats.emp.memory = zeros(size(ERS,1),size(ERS,2),5);
    stats.emp.category = zeros(size(stats.emp.memory));
    stats.perm.memory = zeros(nrep,size(ERS,1),size(ERS,2),5);
    stats.perm.category = zeros(size(stats.perm.memory));
    stats.channels = info.channel_labels;
    
    disp(['computing stats'])
    tic
    for win1 = 1:size(ERS,1) %idx 1 and 2 are now TIME!
        for win2 = 1:size(ERS,2)
            items = squeeze(Reinstatement.ITEM(win1,win2,:));
            
            % memory
            remembered = items(remembered_ix);
            forgotten = items(forgotten_ix);
            [h,p,~,statistics] = ttest2(remembered,forgotten,'alpha',0.05);
            stats.emp.memory(win1,win2,:) = [nanmean(remembered) nanmean(forgotten) h, p, statistics.tstat];
            
            % category
            indoor = items(indoor_ix);
            outdoor = items(outdoor_ix);
            [h,p,~,statistics] = ttest2(indoor,outdoor,'alpha',0.05);
            stats.emp.category(win1,win2,:) = [nanmean(indoor) nanmean(outdoor) h, p, statistics.tstat];
        end
    end
    
   %% PERMUTATION --> Shuffle "item" vector, then repeat all tests from above
    % shuffles item-classification by memory, category etc.. 
    for rep = 1:nrep
        ordered_data = Reinstatement.ITEM;
        shuff_data = ordered_data(:,:,randperm(size(ordered_data,3))); %shuffle once, then compute winXwin with same shuffle!

        for win1 = 1:size(ERS,1) %idx 1 and 2 are now TIME!
            for win2 = 1:size(ERS,2)

                shuffled_items = squeeze(shuff_data(win1,win2,:));

                % memory
                remembered = shuffled_items(remembered_ix);
                forgotten = shuffled_items(forgotten_ix);
                [h,p,~,statistics] = ttest2(remembered,forgotten,'alpha',0.05);
                stats.perm.memory(rep,win1,win2,:) = [nanmean(remembered) nanmean(forgotten) h, p, statistics.tstat];

                % category
                indoor = shuffled_items(indoor_ix);
                outdoor = shuffled_items(outdoor_ix);
                [h,p,~,statistics] = ttest2(indoor,outdoor,'alpha',0.05);
                stats.perm.category(rep,win1,win2,:) = [nanmean(indoor) nanmean(outdoor) h, p, statistics.tstat];
            end
        end
    end
    disp([' took ' num2str(toc) ' seconds'])
    
    saveDir = [dataDir fs fileName fs ersDir fs stats_tag];
    status = mkdir(saveDir);
    
    if save_to_disk
        fileName1 = ['first_level_stats.mat'];
        fileName2 = ['itemspecific_reinstatement.mat'];
        save([saveDir fs fileName1],'stats')
        save([saveDir fs fileName2],'Reinstatement')
    end
end
    
