
function DNN_first_level_btw(subject,dnnType,nrep,corr_type,rsa_tag,roi_tag)
    fs = '\';
    dataDir = '..\DATA';
    ftpath = 'C:\Program Files\MATLAB\R2020b\toolbox\fieldtrip';
        
    %load folder names of subjects
    load('..\subject_list.mat'); %subNames for easier folder accessing
    
    fileName = subNames{subject}; %set folderName of current subject
    disp(fileName)    
    
    %% DNN specifics
    switch dnnType
        case 'AlexNet_avg'
            dnnDir = '.\PlacesNet\AlexNet_Averages';
            layerNum = 8;
            layerNames = {'conv1','conv2','conv3','conv4','conv5','fc6','fc7','fc8'}; %%% ADD FOR OTHERS TOO!
    end
    
    %% load info
    BtwInfo = load([dataDir fs fileName fs 'stats' fs  'BtwItem' fs rsa_tag fs roi_tag fs 'BtwInfo.mat']);
    BtwInfo = BtwInfo.info; 
    numTrials = size(BtwInfo.trialinfo,1);

    %% load DNN pred Mats
    %%%load dnn predictor matrices, (already only half and Fisher Z but not vec -> vectorize
    vec_trial_length = ((numTrials^2) - numTrials) / 2;
    raw_dnns = zeros(layerNum,numTrials,numTrials); %will hold normalized and vectorized lower half dnn matrices

    for layer = 1:layerNum
        if contains(dnnType,'avg')
            load([dataDir fs fileName fs 'stats' fs 'BtwItem' fs 'DNN' fs dnnType fs 'average_' corr_type '_' layerNames{layer} '.mat'])
        end
        raw_dnns(layer,:,:) = RSM;
    end   
    
    %% get indices of trials according to memory / image / complexity etc.
    remembered_ix = find(BtwInfo.trialinfo(:,2) == 1);
    forgotten_ix = find(BtwInfo.trialinfo(:,2) == 2);
    
    %% prepare indices for proper indexing per trial (row-logic but here only half mats!)
    exampleMat = ones(numTrials) - triu(ones(numTrials),1) - eye(numTrials); %half matrix of ones and subj-specific size
    btw_item_indices = zeros(numTrials,numTrials-1); %for each trial (diagonal value), get n Trials - 1 index
    for trial = 1:numTrials
        
        row_idx = find(exampleMat(trial,:));
        col_idx = find(exampleMat(:,trial));
        
        tempMat = exampleMat; tempMat(trial,row_idx) = 10; tempMat(col_idx,trial) = 10;
        linIdx = find(tempMat == 10);
        btw_item_indices(trial,:) = linIdx;
    end
    
    %% load BtwItemSimilarities for the specified phase
    %at this point, all datachans are included
    for phase = 1:2
        tic
        if phase == 1
            disp('ENCODING')
            load([dataDir fs fileName fs 'stats' fs 'BtwItem' fs rsa_tag fs roi_tag fs 'BtwItem_encoding.mat'])
            data = encoding;
        elseif phase == 2
            disp('RETRIEVAL')
            load([dataDir fs fileName fs 'stats' fs 'BtwItem' fs rsa_tag fs roi_tag fs 'BtwItem_recognition.mat'])
            data = recognition;
        end
        
        %%%% STATS
        %%%% 1) Compute item-wise (row) similarity of brain and DNN layer
        %%%% geometries -> yields item-wise similarity score
        %%%% 2) test similarity of remembered vs forgotten per layer
        %%%% 3) repeat test each time shuffling remembered and forgotten!
        %%%% (ERS and BtwItem ERS logic)
        
        
        %% Correlate Brain and DNN data
        
        %get number timepoints
        numWin = size(data,1);
        
        btwitem_dnn_reinstatement = zeros(numTrials,numWin,layerNum);
                
        %% correlate neural data and dnn predictions
        brain_data_flat = reshape(data,numWin,[],1); %make data vector flat (flatten stimxstim since indices are linear!)
        dnn_data_flat = reshape(raw_dnns,layerNum,[],1);
        
        for trial = 1:numTrials
            
            btw_item_trial_brain = brain_data_flat(:,btw_item_indices(trial,:)); %win X btwitemsims (numTrial-1)
            btw_item_trial_dnn = dnn_data_flat(:,btw_item_indices(trial,:)); 
            
            [r1,~] = corr(btw_item_trial_brain',btw_item_trial_dnn','rows','complete','type','spearman'); %should output matrix of size win X Layer
            
            Z1 = 0.5.*log((1+r1)./(1-r1));
            
            btwitem_dnn_reinstatement(trial,:,:) = Z1;
        end
            
        %% now test trial-specific similarities (against 0 and rem vs forg)
        stats = {};
        stats.emp.all = zeros(numWin,layer);
        stats.emp.rem = zeros(numWin,layer);
        stats.emp.forg = zeros(numWin,layer);
        stats.emp.mem = zeros(numWin,layer);
        
        for win = 1:numWin
            for layer = 1:layerNum
                %test reinstatement scores of all against 0
                all_items = squeeze(btwitem_dnn_reinstatement(:,win,layer));
                [~,~,~,stat] = ttest(all_items);
                stats.emp.all(win,layer) = stat.tstat;
                
                %test reinstatment of rem against 0
                rem_items = squeeze(btwitem_dnn_reinstatement(remembered_ix,win,layer));
                [~,~,~,stat] = ttest(rem_items);
                stats.emp.rem(win,layer) = stat.tstat;
                
                %test reinstatment of forg against 0
                forg_items = squeeze(btwitem_dnn_reinstatement(forgotten_ix,win,layer));
                [~,~,~,stat] = ttest(forg_items);
                stats.emp.forg(win,layer) = stat.tstat;
                
                %ttest rem vs forg
                [~,~,~,stat] = ttest2(rem_items,forg_items);
                stats.emp.mem(win,layer) = stat.tstat;
            end
        end
        
        %% permutation ttests, each time shuffling order of trials in btw_item_dnn_reinstatement
        stats.perm.all = zeros(nrep,numWin,layer);
        stats.perm.rem = zeros(nrep,numWin,layer);
        stats.perm.forg = zeros(nrep,numWin,layer);
        stats.perm.mem = zeros(nrep,numWin,layer);
        
        for rep = 1:nrep
            %% now test trial-specific similarities (against 0 and rem vs forg)
            btwitem_dnn_reinstatement_shuffled = btwitem_dnn_reinstatement(randperm(size(btwitem_dnn_reinstatement,1)),:,:);

            for win = 1:numWin
                for layer = 1:layerNum
                    %test reinstatement scores of all against 0
                    all_items = squeeze(btwitem_dnn_reinstatement_shuffled(:,win,layer));
                    [~,~,~,stat] = ttest(all_items);
                    stats.perm.all(rep,win,layer) = stat.tstat;

                    %test reinstatment of rem against 0
                    rem_items = squeeze(btwitem_dnn_reinstatement_shuffled(remembered_ix,win,layer));
                    [~,~,~,stat] = ttest(rem_items);
                    stats.perm.rem(rep,win,layer) = stat.tstat;

                    %test reinstatment of rem against 0
                    forg_items = squeeze(btwitem_dnn_reinstatement_shuffled(forgotten_ix,win,layer));
                    [~,~,~,stat] = ttest(forg_items);
                    stats.perm.forg(rep,win,layer) = stat.tstat;

                    %ttest rem vs forg
                    [~,~,~,stat] = ttest2(rem_items,forg_items);
                    stats.perm.mem(rep,win,layer) = stat.tstat;
                end
            end
        end
        
        %% prep directory and save
        stats.info = BtwInfo;
        
        %set statsdir
        saveDir = [dataDir fs fileName fs 'stats' fs 'BtwItemDNN' fs rsa_tag fs roi_tag];
        status = mkdir(saveDir);

        if phase == 1
            save([saveDir fs 'Item' dnnType '_' corr_type '_enc.mat'],'stats')
            save([saveDir fs 'Item' dnnType '_' corr_type '_enc_rawsim.mat'],'btwitem_dnn_reinstatement')
        elseif phase == 2
            save([saveDir fs 'Item' dnnType '_' corr_type '_rec.mat'],'stats')
            save([saveDir fs 'Item' dnnType '_' corr_type '_rec_rawsim.mat'],'btwitem_dnn_reinstatement')
        end
    end
end


