%% SEA Analysis script
% this is the master script to run the matlab analysis for the SEA Project
% by Nils Kolling, Marius Braunsdorf and Rogier Mars
% This script uses a lot of subfunctions to analyze both behaviour,
% generate Design files for whole brain analyses in FSL and run Timecourse
% analyses (ROI-based, using textfiles compiled with FSL. Both basic Timecourses and PPI's are done)
% To run the script successfully, you will need a couple of additional
% files:
% - 1) subfunctions in the same folder as this master script. also,
% "timecourse_Scripts" (see path specification below)
% - 2) FSL matlab scripts for the timecourse bit of the script [fmt] 
% - 3) behavioural data as log files 
% - 4) folder with timecourses as txt files 
%        Code written by Nils Kolling 2020
%% Set up folder structure
clear all

Data.SubjectIDs = { '01', '02', '05', '06', '07', '08', '09', '10', '11', '12', '13', '14', '16','17','18','19','20','22','23','24','26','27','28','29','32','33','34','35' };
version = '9'; %version of the experimental log files

%% %Define folder locations 
Data.folders.matlaboutput='Matlab_output';
Data.folders.preslog=fullfile('SEAData','Presentation_logfiles');
Data.folders.schedules='scheduleinfo';
Data.Basepath=fullfile('~','Documents','GitHub');
Data.filenames.matlaboutBehave=['behavioural_data_version' version '_' num2str(length(Data.SubjectIDs)) '_subjects'];
Data.folders.timecourses       =fullfile(Data.Basepath,'SEATimecourses','Data');
Data.folders.timecoursesResults=fullfile(Data.Basepath,'SEATimecourses','Results');
%% Read in behavioural data from log files %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DataB.generateDate=date;
DataB.Basepath=Data.Basepath;
for g = 1:length(Data.SubjectIDs) % Subject Loop        
    % Read in files
    readinInfo.filepath=Data.Basepath;
    readinInfo.version = [Data.SubjectIDs{g} '-actual_experiment'];
    readinInfo.filefolder=Data.folders.preslog;
    readinInfo.SubjectIDs=Data.SubjectIDs;
    presinfo(g) = read_presentation_log_sea_mri(readinInfo,g);        
    temp= load(fullfile(Data.Basepath,'SEAData',Data.folders.schedules,[Data.SubjectIDs{g} '_scheduleinfo.mat']));  
    try if isstruct(temp.Schedule)
        Data.schedules(g)=temp.Schedule;
        end
    catch
        if isstruct(temp.scheduleinfo)
            Data.schedules(g)=temp.scheduleinfo;
        else
            keyboard
        end
    end
    
    % Load agent identity in column of matrix (one per Participant)
    DataB=organizeData_makeParameters(presinfo,g,Data.schedules,DataB); % g refers to index in matrices including all participants                
end

%% Add known meta Data  [Is this all the meta data we have and is it correct?]
DataB.meta_data.Subject_Age = [27 25 19 23 33 33 26 24 24 24 19 27 24 21 27 30 21 27 30 25 20 22];
DataB.meta_data.Subject_Gender = {'f', 'f',  'm', 'f','f','m', 'f','f', 'm' ,'f','f','f','f','f','f', 'm', 'm', 'f','f','m','f', 'm'};
DataB.meta_data.Subject_Gender_Binary = [1 1 0 1 1 0 1 1 0 1 1 1 1 1 1 0 0 1 1 0 1 0];

%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Data analysis %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Now run behavioural data analyses, mostly with multiple regressions (glmfit)


%% Self Choice GLM %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Main analysis
DesignMatrixID=1;DataB.GLMs(DesignMatrixID).DesignName='Blue Self Choice using Evidences';    
DataB.GLMs(DesignMatrixID).linkfunction  ='binomial'; %Can be normal, poisson, binomial or gamma
DataB.GLMs(DesignMatrixID).dataTransform ='non';      %Transform can be 'log' or 'non'; 
DataB.GLMs(DesignMatrixID).dataType      ='SelfChoices';      %Could be 'SelfRTs' 'OtherRTs' BothRTs 'SelfCoices' 'OtherChoices' 'BothChoices' 'Ratings'; 
DataB=runRegressions(DataB,DesignMatrixID);
%plotRegressorCorrelations(DataB,DesignMatrixID)   % Optionally by uncommenting this, it is possbile to look at DM regressor correlations
plotGLMResults(DataB,DesignMatrixID,'Effect of regressor on self colour choices (betaweights)')

% % Self and Other choices each separately and sample split
% DesignMatrixID=8;DataB.GLMs(DesignMatrixID).DesignName='Self Choice Samples split';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='binomial';  
% DataB.GLMs(DesignMatrixID).dataTransform ='non';       
% DataB.GLMs(DesignMatrixID).dataType      ='SelfChoices';        
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)
% plotGLMResults(DataB,DesignMatrixID) 

%% Other Choice GLM %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Main analysis
DesignMatrixID=4;DataB.GLMs(DesignMatrixID).DesignName='Blue Other Choice using Evidence';    
DataB.GLMs(DesignMatrixID).linkfunction  ='binomial'; %Can be normal, poisson, binomial or gamma
DataB.GLMs(DesignMatrixID).dataTransform ='non';      %Transform can be 'log' or 'non'; 
DataB.GLMs(DesignMatrixID).dataType      ='OtherChoices';      %Could be 'SelfRTs' 'OtherRTs' BothRTs 'SelfCoices' 'OtherChoices' 'BothChoices' 'Ratings'; 
DataB=runRegressions(DataB,DesignMatrixID);
%plotRegressorCorrelations(DataB,DesignMatrixID)   % Optionally by uncommenting this, it is possbile to look at DM regressor correlations
plotGLMResults(DataB,DesignMatrixID,'Effect of regressor on other colour choices (betaweights)')

% % This analysis shows that occluded bias is only driven by first occluded
% % sample [Don't have this in MS at the moment. Is interesting for occl1 but strange for shared4]
% DesignMatrixID=12;DataB.GLMs(DesignMatrixID).DesignName='Other Choice Samples split';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='binomial';  
% DataB.GLMs(DesignMatrixID).dataTransform ='non';       
% DataB.GLMs(DesignMatrixID).dataType      ='OtherChoices';        
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)
% plotGLMResults(DataB,DesignMatrixID)

% %% Agent specific splits  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %[not in MS at the moment as there was nothing. Howeve, can put in Suppl?]
% %These two analyses show that There is no Agent Specific differences in
% % the use of evidence
% DesignMatrixID=7;DataB.GLMs(DesignMatrixID).DesignName='Other Choice Agent Specific';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='binomial';  
% DataB.GLMs(DesignMatrixID).dataTransform ='non';       
% DataB.GLMs(DesignMatrixID).dataType      ='OtherChoices';        
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)
% plotGLMResults(DataB,DesignMatrixID)  
% 
% DesignMatrixID=10;DataB.GLMs(DesignMatrixID).DesignName='Other Choice Agent x Interaction';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='binomial';  
% DataB.GLMs(DesignMatrixID).dataTransform ='non';       
% DataB.GLMs(DesignMatrixID).dataType      ='OtherChoices';        
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)
% plotGLMResults(DataB,DesignMatrixID) 

% %% Choice Divergence analysis %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% % Here we are using a Choice divergence framework, i.e. is other choice
% % different from self. For sample split it shows a large use of both shared
% % and occluded evidence to determine whether the participant thinks self
% % and other choices should be different
% DesignMatrixID=14;DataB.GLMs(DesignMatrixID).DesignName='Choice Divergence Samples split';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='binomial';  
% DataB.GLMs(DesignMatrixID).dataTransform ='non';       
% DataB.GLMs(DesignMatrixID).dataType      ='ChoiceDivergence';        
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)
% plotGLMResults(DataB,DesignMatrixID)
%  
% DesignMatrixID=15;DataB.GLMs(DesignMatrixID).DesignName='Choice Divergence evid split';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='binomial';  
% DataB.GLMs(DesignMatrixID).dataTransform ='non';       
% DataB.GLMs(DesignMatrixID).dataType      ='ChoiceDivergence';        
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)
% plotGLMResults(DataB,DesignMatrixID)


%% Both Choices GLM %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Main analysis
DesignMatrixID=17;DataB.GLMs(DesignMatrixID).DesignName='Blue Both Choice interaction using Evidences';    
DataB.GLMs(DesignMatrixID).linkfunction  ='binomial'; %Can be normal, poisson, binomial or gamma
DataB.GLMs(DesignMatrixID).dataTransform ='non';      %Transform can be 'log' or 'non'; 
DataB.GLMs(DesignMatrixID).dataType      ='BothChoices';      %Could be 'SelfRTs' 'OtherRTs' BothRTs 'SelfCoices' 'OtherChoices' 'BothChoices' 'Ratings'; 
DataB=runRegressions(DataB,DesignMatrixID);
%plotRegressorCorrelations(DataB,DesignMatrixID)
plotGLMResults(DataB,DesignMatrixID,'Effect of regressor on self & other colour choices (betaweights)')

DesignMatrixID=18;DataB.GLMs(DesignMatrixID).DesignName='Blue Both Choice separate regressors using Evidences';    
DataB.GLMs(DesignMatrixID).linkfunction  ='binomial'; %Can be normal, poisson, binomial or gamma
DataB.GLMs(DesignMatrixID).dataTransform ='non';      %Transform can be 'log' or 'non'; 
DataB.GLMs(DesignMatrixID).dataType      ='BothChoices';      %Could be 'SelfRTs' 'OtherRTs' BothRTs 'SelfCoices' 'OtherChoices' 'BothChoices' 'Ratings'; 
DataB=runRegressions(DataB,DesignMatrixID);
%plotRegressorCorrelations(DataB,DesignMatrixID)
plotGLMResults(DataB,DesignMatrixID,'Effect of regressor on self & other colour choices (betaweights)')


% %% Supplemental analyses [not in MS]
% DesignMatrixID=5;DataB.GLMs(DesignMatrixID).DesignName='Other Absolute Evidences RTs';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='gamma'; 
% DataB.GLMs(DesignMatrixID).dataTransform ='log';      
% DataB.GLMs(DesignMatrixID).dataType      ='OtherRTs';       
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)
% plotGLMResults(DataB,DesignMatrixID)
% 
% DesignMatrixID=6;DataB.GLMs(DesignMatrixID).DesignName='Other Accuracy using evidences';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='binomial';  
% DataB.GLMs(DesignMatrixID).dataTransform ='non';       
% DataB.GLMs(DesignMatrixID).dataType      ='OtherAccuracy';        
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)
% plotGLMResults(DataB,DesignMatrixID)  
% 
% % Supplemental Analyses of RT and choices in terms of accuracy [not reported in MS]  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% DesignMatrixID=2;DataB.GLMs(DesignMatrixID).DesignName='Self Absolute Evidences RTs';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='gamma'; 
% DataB.GLMs(DesignMatrixID).dataTransform ='log';      
% DataB.GLMs(DesignMatrixID).dataType      ='SelfRTs';       
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)   % Optionally by uncommenting this, it is possbile to look at DM regressor correlations
% plotGLMResults(DataB,DesignMatrixID)
% 
% DesignMatrixID=3;DataB.GLMs(DesignMatrixID).DesignName='Self Accuracy using evidences';    
% DataB.GLMs(DesignMatrixID).linkfunction  ='binomial';  
% DataB.GLMs(DesignMatrixID).dataTransform ='non';       
% DataB.GLMs(DesignMatrixID).dataType      ='SelfAccuracy';        
% DataB=runRegressions(DataB,DesignMatrixID);
% %plotRegressorCorrelations(DataB,DesignMatrixID)   % Optionally by uncommenting this, it is possbile to look at DM regressor correlations
% plotGLMResults(DataB,DesignMatrixID)  

%% Other plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
plotBehaviours(DataB)

%% generating FSL glm Design matrices here %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fMRIGLM_info.fMRIDesignName      ='Sampling_ChosenEvid';
fMRIGLM_info.SelfChoiceSelection ='AllChoicesSelf';
fMRIGLM_info.OtherChoiceSelection='AllChoicesOther';
fMRIGLM_info.SamplingSelection   ='AllSamples';
fMRIGLM_info.Outcomeselection    ='AllOutcomes';
[DataB,Sampling_ChosenEvid]=generate_fMRI_GLMS(DataB,fMRIGLM_info);

%% Timecourse analysis scripts %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%These scripts run the timecourse analysis based on ROI timecourses
%extracted using FSL (both basic and PPI's controlling for global timecourse in the regression as additional regressor)
addpath(fullfile(Data.Basepath, 'git_repository', 'timecourse_Scripts'));% Main analysis for MS
[Data,Results]=RunTimeCourse(Sampling_ChosenEvid,Data,'Sampling_ChosenEvid',DataB);



%% Only other additional analysis, potentially worth discussing, but currently not included in manuscript. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% fMRIGLM_info.fMRIDesignName      ='ChoiceDivergenceAnalysis';
% fMRIGLM_info.SelfChoiceSelection ='AllChoicesSelf';
% fMRIGLM_info.OtherChoiceSelection='AllChoicesOther';
% fMRIGLM_info.SamplingSelection   ='AllSamples';
% fMRIGLM_info.Outcomeselection    ='AllOutcomes';
% [DataB,ChoiceDivergenceGLM]=generate_fMRI_GLM_extensive(DataB,fMRIGLM_info);
% 
% % there is some stuff in here about occluded counter evidence in the
% % sampling period and something about absEvidence, particularly occluded
% % counterevidence being in Fpole. occl counterevidence is different from
% % other analysis "Sampling_ChosenEvid" because evidence is not coded for
% % choice. Interestingly without choice coding all the dmPFC stuff
% % disappears. Self-Choice however looks identical to other framing i.e.
% % directions are as one would predict from choice analysis.
% [Data,Results]=RunTimeCourse(ChoiceDivergenceGLM,Data,'ChoiceDivergenceAnalysis',DataB);
% 

keyboard
 

