add stateswitch preprocessing scripts

code/01_prepare_preprocessing.m

@@ -0,0 +1,189 @@
+%% A_prepare_preprocessing_StateSwitch_170913
+
+% prepare (i.e. filter + downsample data for ICA1)
+
+% INPUT: eeglab structures: EEG data, merged with eye tracking data
+
+% 170913 | JQK adapted function from MD's LC-NE scripts
+% 180205 | added OA IDs, corrected paths
+% 180207 | 2227 changed channels
+
+%% initialize
+
+restoredefaultpath;
+clear all; close all; pack; clc;
+
+%% pathdef
+
+%pn.eeg_root     = '/Users/kosciessa/Desktop/mountpoint_tardis_LNDG/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/';
+pn.study        = '/Volumes/LNDG/Projects/StateSwitch/';
+pn.eeg_root     = [pn.study, 'dynamic/data/eeg/task/A_preproc/SA_preproc_study/'];
+pn.dynamic_In   = [pn.eeg_root, 'B_data/B_EEG_ET_ByRun/'];
+pn.eyePlots     = [pn.eeg_root, 'C_figures/B_EyeChannels/']; mkdir(pn.eyePlots);
+pn.EEG_out      = [pn.eeg_root, 'B_data/C_EEG_FT/']; mkdir(pn.EEG_out);
+pn.History_out  = [pn.eeg_root, 'B_data/D_History/']; mkdir(pn.History_out);
+        
+%% add ConMemEEG tools
+
+pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB'];            addpath(genpath(pn.MWBtools));
+pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG'];            addpath(genpath(pn.THGtools));
+pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common'];         addpath(genpath(pn.commontools));
+pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+pn.scripts      = [pn.eeg_root, 'A_scripts/helper/']; addpath(pn.scripts);
+
+%% datadef
+
+condEEG = 'dynamic';
+
+%% define IDs for preprocessing
+
+% % N = 47 (1213 excluded);
+% IDs = {'1117';'1118';'1120';'1124';'1126';'1131';'1132';'1135';'1136';'1138';'1144';'1151';'1158';'1160';'1163';'1164';'1167';'1169';'1172';'1173';'1178';'1182';'1215';'1216';'1219';'1221';'1223';'1227';'1228';'1233';'1234';'1237';'1239';'1240';'1243';'1245';'1247';'1250';'1252';'1257';'1261';'1265';'1266';'1268';'1270';'1276';'1281'};
+
+% N = 53 OAs;
+IDs = {'2104';'2107';'2108';'2112';'2118';'2120';'2121';'2123';'2125';'2129';'2130';'2131';'2132';'2133';'2134';'2135';'2139';'2140';'2145';'2147';'2149';'2157';'2160';'2201';'2202';'2203';'2205';'2206';'2209';'2210';'2211';'2213';'2214';'2215';'2216';'2217';'2219';'2222';'2224';'2226';'2227';'2236';'2237';'2238';'2241';'2244';'2246';'2248';'2250';'2251';'2252';'2258';'2261'};
+
+%% loop IDs
+for id = 1:length(IDs)
+
+%%  preprocessing
+%if ~exist([pdat.analyses.EEG_dat '\' num2str(ID(id)) '\raw_eeg\' num2str(ID(id)) '_' condEEG '_EEG_Raw_Rlm_Flh_Res.mat'],'file')
+
+%     try
+    
+    for indRun = 1:4
+    
+        display(['processing ID ' IDs{id}, ' Run ', num2str(indRun)]);
+        
+        file = [pn.([condEEG, '_In']), IDs{id} ,'_r',num2str(indRun), '_', condEEG, '_eyeEEG.set'];
+
+    %%  load data
+
+        cfg_preprocessing           = [];
+        cfg_preprocessing.datafile  = file;
+
+        % load header & event information
+        config              = ft_read_header(cfg_preprocessing.datafile);
+        config.data_file    = cfg_preprocessing.datafile;
+        config.mrk          = ft_read_event(cfg_preprocessing.datafile);
+
+        % define reading & preprocessing parameters
+        cfg_preprocessing.channel     = {'all'};
+        cfg_preprocessing.implicitref = 'REF'; % recover implicit reference
+        
+        % get all data first, then apply specific steps only for subsets
+        data_eyeEEG = ft_preprocessing(cfg_preprocessing);
+
+        %% SWITCH CHANNELS ACCORDING TO ARRANGEMENT!
+        
+%         if max(strcmp(IDs{id}, {'1213'}))==1
+%             data_eyeEEG = SS_switchChannels_Pilot_Study(data_eyeEEG);
+        if max(strcmp(IDs{id}, {'1126'; '2227'}))==1
+            data_eyeEEG = SS_switchChannels_GreenYellow(data_eyeEEG);       % green and yellow boxes exchanged
+            data_eyeEEG = SS_switchChannels_Study_noA1(data_eyeEEG);        % TP9 and TP10 exchanged manually
+        elseif max(strcmp(IDs{id}, {'1216'}))==1
+            data_eyeEEG = SS_switchChannels_GreenYellow(data_eyeEEG);       % green and yellow boxes exchanged
+            data_eyeEEG = SS_switchChannels_Study(data_eyeEEG);             % TP9(A1) and TP10(FCz) wrongly ordered in workspace
+        elseif max(strcmp(IDs{id}, {'1118'; '1215'; '1124'}))==1
+            data_eyeEEG = SS_switchChannels_Study_noA1(data_eyeEEG);        % TP9 and TP10 exchanged manually
+        else 
+            data_eyeEEG = SS_switchChannels_Study(data_eyeEEG);             % TP9(A1) and TP10(FCz) wrongly ordered in workspace
+        end
+        data_eyeEEG = SS_changeChannels(data_eyeEEG);
+        
+        %% Plot excerpt from eye channels
+        
+        h = figure;
+        subplot(3,1,1);
+        plot(zscore(data_eyeEEG.trial{1,1}(63,5000:15000),[],2)); % LHEOG
+        hold on; plot(zscore(data_eyeEEG.trial{1,1}(64,5000:15000),[],2)); % RHEOG
+        xlim([0 10000]); legend({'LHEOG'; 'RHEOG'});
+        subplot(3,1,2);
+        %hold on; plot(zscore(data_eyeEEG.trial{1,1}(65,5000:15000),[],2)); % ECG
+        hold on; plot(zscore(data_eyeEEG.trial{1,1}(62,5000:15000),[],2)); % IOR
+        hold on; plot(zscore(data_eyeEEG.trial{1,1}(24,5000:15000),[],2)); % FT8
+        xlim([0 10000]); legend({'IOR'; 'FT8'});
+        subplot(3,1,3);
+        if size(data_eyeEEG.trial{1,1},1)>66
+            hold on; plot(zscore(data_eyeEEG.trial{1,1}(67,5000:15000),[],2)); % left Gaze X
+            hold on; plot(zscore(data_eyeEEG.trial{1,1}(68,5000:15000),[],2)); % left Gaze Y
+            xlim([0 10000]);legend({'left X'; 'left Y'});
+        end
+        pn.plotFolder = pn.eyePlots;
+        figureName = [IDs{id} ,'_r',num2str(indRun), '_', condEEG, '_eyeChannels'];
+        saveas(h, [pn.plotFolder, figureName], 'png');
+        close(h);
+
+        %% continue preprocessing
+        
+        % spare ECG from reref!
+        tmp.data_ECG = data_eyeEEG;
+
+        cfg_ECG                 = [];
+        cfg_ECG.channel         = {'ECG'};
+        cfg_ECG.bsfilter        = 'yes';
+        cfg_ECG.bsfiltord       = 4;
+        cfg_ECG.bsfreq          = [48 52];
+        tmp.data_ECG            = ft_preprocessing(cfg_ECG, tmp.data_ECG);
+
+        % exclude ET, ECG channel (and others) here?
+
+        cfg_preprocessing               = [];                                 
+        cfg_preprocessing.channel       = {'all','-TIME','-L-GAZE-X','-L-GAZE-Y','-L-AREA','-L-VEL-X','-L-VEL-Y','-RES-X','-RES-Y','-INPUT'};
+
+        cfg_preprocessing.continuous  = 'yes';
+        cfg_preprocessing.demean      = 'yes';
+
+        cfg_preprocessing.reref       = 'yes';
+        cfg_preprocessing.refchannel  = {'A1', 'A2'};
+        cfg_preprocessing.implicitref = 'A2';
+
+        cfg_preprocessing.hpfilter    = 'yes';
+        cfg_preprocessing.hpfreq      = 1;
+        cfg_preprocessing.hpfiltord   = 4;
+        cfg_preprocessing.hpfilttype  = 'but';
+
+        cfg_preprocessing.lpfilter    = 'yes';
+        cfg_preprocessing.lpfreq      = 100;
+        cfg_preprocessing.lpfiltord   = 4;
+        cfg_preprocessing.lpfilttype  = 'but';
+
+    %     cfg_preprocessing.bsfilter        = 'yes';
+    %     cfg_preprocessing.bsfiltord       = 4;
+    %     cfg_preprocessing.bsfreq          = [48 52];
+
+        % define settings for resampling
+        cfg_resample.resamplefs = 250;
+        cfg_resample.detrend    = 'no';
+        cfg_resample.feedback   = 'no';
+        cfg_resample.trials     = 'all';
+
+        % get data
+        % load complete data at once & resample afterwards
+        data_EEG = ft_preprocessing(cfg_preprocessing, data_eyeEEG);
+
+        % after reref, copy original ECG (we don't want it rerefed)
+        data_EEG.trial{1,1}(find(strcmp(data_EEG.label , 'ECG')),:) = tmp.data_ECG.trial{1,1};
+        clear tmp;
+
+        data_EEG = ft_resampledata(cfg_resample,data_EEG);
+        % change data precision to single
+        for t = 1:length(data_EEG.trial)
+            data_EEG.trial{t} = single(data_EEG.trial{t});
+        end; clear t
+
+    %%  remove data at edges
+        % not implemented & maybe not necessary
+    %%  save data
+
+        save([pn.EEG_out, IDs{id}, '_r',num2str(indRun), '_', condEEG, '_EEG_Raw_Rlm_Flh_Res.mat'],'data_EEG', '-v7.3');
+        save([pn.EEG_out, IDs{id}, '_r',num2str(indRun), '_', condEEG, '_eyeEEG_Raw.mat'],'data_eyeEEG', '-v7.3');
+        % save config
+        save([pn.History_out, IDs{id}, '_r',num2str(indRun), '_', condEEG, '_config.mat'],'config');
+        % clear variables
+        clear cfg_* config data_EEG data_eyeEEG file
+    
+    end % run loop
+
+end; clear id;

code/02_visual_inspection.m

@@ -0,0 +1,117 @@
+%% B_visual_inspection_for_ica_170913
+
+% 170913 | JQK adapted from MD script
+% 180207 | adapted for STSWD OAs
+
+%% initialize
+
+restoredefaultpath;
+clear all; close all; pack; clc;
+
+%% pathdef
+
+pn.study        = '/Volumes/LNDG/Projects/StateSwitch/';
+pn.eeg_root     = [pn.study, 'dynamic/data/eeg/task/A_preproc/SA_preproc_study/'];
+%pn.eeg_root = '/Users/kosciessa/Desktop/mountpoint_tardis_LNDG/LNDG/StateSwitch/WIP_eeg/SA_preproc_study/';
+pn.dynamic_In   = [pn.eeg_root, 'B_data/B_EEG_ET_ByRun/'];
+pn.eyePlots     = [pn.eeg_root, 'C_figures/B_EyeChannels/']; mkdir(pn.eyePlots);
+pn.EEG          = [pn.eeg_root, 'B_data/C_EEG_FT/']; mkdir(pn.EEG);
+pn.History      = [pn.eeg_root, 'B_data/D_History/']; mkdir(pn.History);
+
+%% add ConMemEEG tools
+
+pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB'];            addpath(genpath(pn.MWBtools));
+pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG'];            addpath(genpath(pn.THGtools));
+pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common'];         addpath(genpath(pn.commontools));
+pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+
+%% define Condition & IDs for preprocessing
+
+condEEG = 'dynamic';
+
+%% define IDs for visual screening
+
+% N = 48;
+IDs = {'1117';'1118';'1120';'1124';'1126';'1131';'1132';'1135';'1136';'1138';'1144';'1151';'1158';'1160';'1163';'1164';'1167';'1169';'1172';'1173';'1178';'1182';'1213';'1215';'1216';'1219';'1221';'1223';'1227';'1228';'1233';'1234';'1237';'1239';'1240';'1243';'1245';'1247';'1250';'1252';'1257';'1261';'1265';'1266';'1268';'1270';'1276';'1281'};
+
+% N = 53 OAs;
+IDs = {'2104';'2107';'2108';'2112';'2118';'2120';'2121';'2123';'2125';'2129';'2130';'2131';'2132';'2133';'2134';'2135';'2139';'2140';'2145';'2147';'2149';'2157';'2160';'2201';'2202';'2203';'2205';'2206';'2209';'2210';'2211';'2213';'2214';'2215';'2216';'2217';'2219';'2222';'2224';'2226';'2227';'2236';'2237';'2238';'2241';'2244';'2246';'2248';'2250';'2251';'2252';'2258';'2261'};
+
+%% loop IDs
+
+for id = 4:length(IDs)
+    for indRun = 1:4
+        %% load data, start screening
+        clc;
+        ID = str2num(IDs{id}); disp(['Processing ',IDs{id}, ' Run ', num2str(indRun)]);
+        % load config
+        load([pn.History, num2str(ID), '_r',num2str(indRun), '_', condEEG, '_config.mat'],'config')
+        % load data
+        load([pn.EEG, num2str(ID), '_r',num2str(indRun), '_', condEEG, '_EEG_Raw_Rlm_Flh_Res.mat'],'data_EEG')
+        % data browser
+        if ~isfield(config,'visual_inspection')
+            cfg = [];
+            cfg.continuous      = 'yes';                           
+            cfg.preproc.demean  = 'yes';  % Baselinekorrektur
+            cfg.blocksize       = 125;
+            cfg.ylim            = [-50 50];
+            cfg.viewmode        = 'vertical';
+            cfg.plotlabels      = 'yes';
+            % cfg.channelcolormap = [];
+        else
+            cfg = config.visual_inspection;
+        end
+
+        cfg.channel = {'all','-ECG', '-A1', '-REF'};
+
+        % inspect data
+        cfg = ft_databrowser(cfg,data_EEG);
+        
+        % check if ICA conducted on this dataset; if so new artifacts not saved
+        if ~isfield(config,'ica1')
+
+            % marked segments
+            reject = cfg.artfctdef.visual.artifact;
+
+            trl = [];
+            
+            % generate "trial" (i.e. segment) structure for ICA
+            if exist('reject','var')
+                if min(size(reject)) ~= 0 % JQK fix: 170915
+                    for j = 1:size(reject,1)+1
+                        if j == 1
+                            trl(j,:) = [1                reject(j,1)-1  0];
+                        elseif j <= size(reject,1)
+                            trl(j,:) = [reject(j-1,2)+1  reject(j,1)-1  0];
+                        else
+                            trl(j,:) = [reject(j-1,2)+1  size(data_EEG.trial{1},2) 0];
+                        end
+                    end; clear j
+                else
+                    trl = [1 size(data_EEG.trial{1},2) 0];
+                end
+            end
+
+            % eliminate trials without data points
+            ex = find((trl(:,1)-trl(:,2))==0);
+            trl(ex,:) = []; clear ex
+
+            % eliminate trials shorter than 3 sec (! srate = 250 Hz)
+            ex = find((trl(:,2)-trl(:,1))<750);
+            trl(ex,:) = []; clear ex
+
+            % add to config structure
+            config.visual_inspection = cfg;
+            config.trl_ica1          = trl;
+
+            % save config
+            save([pn.History, num2str(ID), '_r',num2str(indRun), '_', condEEG, '_config.mat'],'config')
+
+        else
+
+            warning('ICA already conducted. Changes not saved.')
+
+        end
+    end % run loop
+end % id loop

code/03_ica/E_STSW_ica1_180108.m

@@ -0,0 +1,175 @@
+function E_STSW_ica1_180108(id)
+
+% 170915 | JQK adapted from MD script
+% 180103 | adapted for STSW Study
+% 180108 | adapted for tardis
+
+%% initialize
+
+%restoredefaultpath;
+%clear all; close all; pack; clc;
+
+%% pathdef
+
+if ismac
+    pn.study        = '/Volumes/LNDG/Projects/StateSwitch/';
+    pn.eeg_root     = [pn.study, 'dynamic/data/eeg/task/A_preproc/SA_preproc_study_YA/'];
+    pn.dynamic_In   = [pn.eeg_root, 'B_data/B_EEG_ET_ByRun/'];
+    pn.eyePlots     = [pn.eeg_root, 'C_figures/B_EyeChannels/']; mkdir(pn.eyePlots);
+    pn.EEG          = [pn.eeg_root, 'B_data/C_EEG_FT/']; mkdir(pn.EEG);
+    pn.History      = [pn.eeg_root, 'B_data/D_History/']; mkdir(pn.History);
+    % add ConMemEEG tools
+    pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB/'];            addpath(genpath(pn.MWBtools));
+    pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG/'];            addpath(genpath(pn.THGtools));
+    pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common/'];         addpath(genpath(pn.commontools));
+    pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+    pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+else
+    pn.root         = '/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/';
+    pn.EEG          = [pn.root, 'B_data/C_EEG_FT/']; mkdir(pn.EEG);
+    pn.History      = [pn.root, 'B_data/D_History/']; mkdir(pn.History);
+    pn.THGtools     = [pn.root, 'T_tools/fnct_THG/'];
+    pn.commontools  = [pn.root, 'T_tools/fnct_common/'];
+    % add external tools (need to be compiled in)
+    %pn.commontools  = [pn.root, 'T_Tools/fnct_common'];         addpath(genpath(pn.commontools));
+    %pn.FT           = [pn.root, 'T_Tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+end
+
+%% define Condition & IDs for preprocessing
+
+condEEG = 'dynamic';
+
+%% define IDs for visual screening
+
+% N = 48;
+IDs = {'1117';'1118';'1120';'1124';'1126';'1131';'1132';'1135';'1136';'1138';'1144';'1151';'1158';'1160';'1163';'1164';'1167';'1169';'1172';'1173';'1178';'1182';'1213';'1215';'1216';'1219';'1221';'1223';'1227';'1228';'1233';'1234';'1237';'1239';'1240';'1243';'1245';'1247';'1250';'1252';'1257';'1261';'1265';'1266';'1268';'1270';'1276';'1281'};
+
+id = str2num(id);
+
+%% loop IDs
+%for id = 1:length(IDs)
+	display(['processing ID ' num2str(IDs{id})]);
+    try
+        if 1%~exist([pn.EEG, IDs{id}, '_r',num2str(indRun), '_', condEEG, '_EEG_Rlm_Fhl_Ica.mat'],'file')
+            
+            for iRun = 1:4
+                %%  load raw data & exclude parts containing artifacts
+
+                % load config
+                config = [];
+                load([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config.mat'],'config');
+
+                config_tmp.(['trl_', num2str(iRun)]) = config.trl_ica1;
+                config_tmp.(['visual_inspection_', num2str(iRun)]) = config.visual_inspection;
+
+                % define segment(s) to be read by fieldtrip
+                cfg.trl = config.trl_ica1;
+
+                % load data
+                dataByRun{iRun} = load([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_EEG_Raw_Rlm_Flh_Res.mat'],'data_EEG');
+
+                % adjust final length
+                cfg.trl(end,2) = size(dataByRun{iRun}.data_EEG.trial{1},2);
+                
+                % "segment" data
+                dataByRun{iRun}.data = cm_segmentation_of_continuous_data_fieldtrip_20150825(dataByRun{iRun}.data_EEG,cfg);    
+
+                % clear cfg structure
+                clear cfg
+
+                %%  segmentation for ICA
+
+                % define settings
+                cfg.length = 2;
+                cfg.n      = 5000;      % keep all possible trials
+                cfg.type   = 'rnd';     % select trials randomly if > 1500 trials available
+                cfg.seed   = 20170915 + str2num(IDs{id});
+
+                % arbitrary segmentation - segments a 2 sec
+                % NOTE original segments will be overwritten
+                dataByRun{iRun}.data = cm_arbitrary_segmentation_fieldtrip_20150210(dataByRun{iRun}.data,cfg);
+                % data = ft_checkdata(data,'feedback','yes');
+            end
+
+            % combine data into one structure for ICA
+            data = ft_appenddata(cfg, dataByRun{1}.data, dataByRun{2}.data, dataByRun{3}.data, dataByRun{4}.data);
+
+            %%  ICA
+
+            % date
+            dt = date;
+
+            % ica config
+            cfg.method           = 'runica';
+            cfg.channel          = {'all','-ECG','-A2'};                % additional channel should be excluded already...
+            cfg.trials           = 'all';
+            cfg.numcomponent     = 'all';
+            cfg.demean           = 'no';
+            cfg.runica.extended  = 1;
+            cfg.runica.logfile   = [pn.History 'log_' IDs{id} '_' condEEG '_ICA1_' dt '.txt'];
+
+            % run ICA
+            icadat = ft_componentanalysis(cfg,data);
+
+            %%  automatic ICA labeling
+
+            [iclabels] = cm_automatic_IC_detection_20170919(data,icadat);
+
+            %%  save data for ICA labeling
+
+            config.trl_1 = config_tmp.trl_1;
+            config.trl_2 = config_tmp.trl_2;
+            config.trl_3 = config_tmp.trl_3;
+            config.trl_4 = config_tmp.trl_4;
+
+            config.visual_inspection_1 = config_tmp.visual_inspection_1;
+            config.visual_inspection_2 = config_tmp.visual_inspection_2;
+            config.visual_inspection_3 = config_tmp.visual_inspection_3;
+            config.visual_inspection_4 = config_tmp.visual_inspection_4;
+
+            % - include ICA solution in data
+            data.topo 	   = icadat.topo;
+            data.unmixing  = icadat.unmixing;
+            data.topolabel = icadat.topolabel;
+            data.cfg       = icadat.cfg;
+
+            % - include ICA solution in config
+            config.ica1.date      = dt;
+            config.ica1.topo 	  = icadat.topo;
+            config.ica1.unmixing  = icadat.unmixing;
+            config.ica1.topolabel = icadat.topolabel;
+            config.ica1.cfg       = icadat.cfg;
+            config.ica1.iclabels.auto = iclabels;
+
+            % fieldtrip format electrode information
+            load([pn.THGtools, 'electrodelayouts/realistic_1005.mat'])
+            data.elec = cm_elec2dataelec_20170919(realistic_1005,data);
+
+            % EEGLAB format electrode information
+            load([pn.THGtools 'electrodelayouts/chanlocs_eeglab_MPIB_64_electrodes.mat'])
+            data = cm_chanlocs2MPIB64_20140126(data,chanlocs);
+
+            % - include channel information in config
+            config.elec     = data.elec;
+            config.chanlocs = data.chanlocs;
+
+            % keep ICA labels
+            data.iclabels = iclabels;
+
+            % save data
+            save([pn.EEG, IDs{id}, '_', condEEG, '_EEG_Rlm_Fhl_Ica.mat'],'data')
+
+            % save config
+            config.preproc_version = '20170915';
+            save([pn.History, IDs{id}, '_', condEEG, '_config.mat'],'config')
+
+            % clear variables
+            clear chanlocs cfg config data dt icadat iclabels realistic_1005
+
+         end % file available
+    catch ME
+        warning('Error occured. Please check.');
+        rethrow(ME)
+    end % try..catch
+%end; clear id
+

code/03_ica/E_STSW_ica1_180108_START.sh

@@ -0,0 +1,22 @@
+#!/bin/bash
+
+# call the BOSC analysis by session and participant
+cd /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/E_STSW_ica1_180108/
+
+IDs=(1117 1118 1120 1124 1125 1126 1131 1132 1135 1136 1151 1160 1164 1167 1169 1172 1173 1178 1182 1214 1215 1216 1219 1223 1227 1228 1233 1234 1237 1239 1240 1243 1245 1247 1250 1252 1257 1261 1265 1266 1268 1270 1276 1281 2142 2253 2254 2255)
+
+mkdir /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/E_STSW_ica1_180108
+
+for i in $(seq 1 ${#IDs[@]}); do
+	echo "#PBS -N STSW_eeg_ica_${i}" 										>> job
+	echo "#PBS -l walltime=03:00:00" 										>> job
+	echo "#PBS -l mem=8gb" 													>> job
+	echo "#PBS -j oe" 														>> job
+	echo "#PBS -o /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/E_STSW_ica1_180108" >> job
+	echo "#PBS -m n" 														>> job
+	echo "#PBS -d ." 														>> job
+	echo "./E_STSW_ica1_180108_run.sh /opt/matlab/R2016b $i " 			>> job
+	qsub job
+	rm job
+done
+done

code/03_ica/E_STSW_ica1_180108_START_repeat.sh

@@ -0,0 +1,31 @@
+#!/bin/bash
+
+# call the BOSC analysis by session and participant
+cd /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/E_STSW_ica1_180108/
+
+IDs=(1117 1118 1120 1124 1125 1126 1131 1132 1135 1136 1151 1160 1164 1167 1169 1172 1173 1178 1182 1214 1215 1216 1219 1223 1227 1228 1233 1234 1237 1239 1240 1243 1245 1247 1250 1252 1257 1261 1265 1266 1268 1270 1276 1281 2142 2253 2254 2255)
+
+mkdir /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/E_STSW_ica1_180108
+
+echo "#PBS -N STSW_eeg_ica_7" 											>> job
+echo "#PBS -l walltime=30:00:00" 										>> job
+echo "#PBS -l mem=8gb" 													>> job
+echo "#PBS -j oe" 														>> job
+echo "#PBS -o /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/E_STSW_ica1_180108" >> job
+echo "#PBS -m n" 														>> job
+echo "#PBS -d ." 														>> job
+echo "./E_STSW_ica1_180108_run.sh /opt/matlab/R2016b 7 " 				>> job
+qsub job
+rm job
+
+
+echo "#PBS -N STSW_eeg_ica_34" 											>> job
+echo "#PBS -l walltime=30:00:00" 										>> job
+echo "#PBS -l mem=8gb" 													>> job
+echo "#PBS -j oe" 														>> job
+echo "#PBS -o /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/E_STSW_ica1_180108" >> job
+echo "#PBS -m n" 														>> job
+echo "#PBS -d ." 														>> job
+echo "./E_STSW_ica1_180108_run.sh /opt/matlab/R2016b 34 " 				>> job
+qsub job
+rm job

code/03_ica/E_STSW_ica1_180108_START_repeat2.sh

@@ -0,0 +1,19 @@
+#!/bin/bash
+
+# call the BOSC analysis by session and participant
+cd /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/E_STSW_ica1_180108/
+
+IDs=(1117 1118 1120 1124 1125 1126 1131 1132 1135 1136 1151 1160 1164 1167 1169 1172 1173 1178 1182 1214 1215 1216 1219 1223 1227 1228 1233 1234 1237 1239 1240 1243 1245 1247 1250 1252 1257 1261 1265 1266 1268 1270 1276 1281 2142 2253 2254 2255)
+
+mkdir /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/E_STSW_ica1_180108
+
+echo "#PBS -N STSW_eeg_ica_21" 											>> job
+echo "#PBS -l walltime=15:00:00" 										>> job
+echo "#PBS -l mem=8gb" 													>> job
+echo "#PBS -j oe" 														>> job
+echo "#PBS -o /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/E_STSW_ica1_180108" >> job
+echo "#PBS -m n" 														>> job
+echo "#PBS -d ." 														>> job
+echo "./E_STSW_ica1_180108_run.sh /opt/matlab/R2016b 21 " 				>> job
+qsub job
+rm job

code/03_ica/E_STSW_ica1_180108_prepare.sh

@@ -0,0 +1,22 @@
+#!/bin/bash
+
+# This script prepares tardis by compiling the necessary function in MATLAB.
+
+#ssh tardis # access tardis
+
+# check and choose matlab version
+#module avail matlab
+module load matlab/R2016b
+
+# compile functions
+
+matlab
+%% add fieldtrip toolbox
+addpath('/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/T_tools/fieldtrip-20170904/')
+ft_defaults()
+ft_compile_mex(true)
+%% go to analysis directory containing .m-file
+cd('/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/E_STSW_ica1_180108/')
+%% compile function and append dependencies
+mcc -m E_STSW_ica1_180108.m -a /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/T_tools/fnct_common -a /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/T_tools/eeglab14_1_1b/functions/sigprocfunc
+exit

code/03_ica/E_STSW_ica1_180108_run.sh

@@ -0,0 +1,33 @@
+#!/bin/sh
+# script for execution of deployed applications
+#
+# Sets up the MATLAB Runtime environment for the current $ARCH and executes 
+# the specified command.
+#
+exe_name=$0
+exe_dir=`dirname "$0"`
+echo "------------------------------------------"
+if [ "x$1" = "x" ]; then
+  echo Usage:
+  echo    $0 \<deployedMCRroot\> args
+else
+  echo Setting up environment variables
+  MCRROOT="$1"
+  echo ---
+  LD_LIBRARY_PATH=.:${MCRROOT}/runtime/glnxa64 ;
+  LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/bin/glnxa64 ;
+  LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/sys/os/glnxa64;
+  LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/sys/opengl/lib/glnxa64;
+  export LD_LIBRARY_PATH;
+  echo LD_LIBRARY_PATH is ${LD_LIBRARY_PATH};
+  shift 1
+  args=
+  while [ $# -gt 0 ]; do
+      token=$1
+      args="${args} \"${token}\"" 
+      shift
+  done
+  eval "\"${exe_dir}/E_STSW_ica1_180108\"" $args
+fi
+exit
+

code/03_ica/mccExcludedFiles.log

@@ -0,0 +1,10 @@
+The List of Excluded Files
+Excluded files	 Exclusion Message ID	 Reason For Exclusion	 Exclusion Rule
+/opt/matlab/R2016b/toolbox/local/+matlab/+internal/+toolboxes/addInstalledToolboxesToJavaPath.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/[+]matlab/[+]internal/[+]toolboxes/
+/opt/matlab/R2016b/toolbox/local/+matlab/+internal/+toolboxes/addInstalledToolboxesToPath.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/[+]matlab/[+]internal/[+]toolboxes/
+/opt/matlab/R2016b/toolbox/local/+matlab/+internal/+zipAddOns/addInstalledZipAddOnsToPath.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/[+]matlab/[+]internal/[+]zipAddOns/
+/opt/matlab/R2016b/toolbox/local/pathdef.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/pathdef[.]m
+/opt/matlab/R2016b/toolbox/matlab/codetools/@mtree/Arg.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/matlab/codetools/
+/opt/matlab/R2016b/toolbox/matlab/codetools/@mtree/Index.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/matlab/codetools/
+/opt/matlab/R2016b/toolbox/matlab/codetools/@mtree/mtree.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/matlab/codetools/
+/opt/matlab/R2016b/toolbox/matlab/codetools/initdesktoputils.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/matlab/codetools/

code/04_ica_labeling.m

@@ -0,0 +1,99 @@
+%% D_ICA_labeling_resting_20141217
+
+% The central idea here is to take the full data and label the ICA based on
+% the data across runs. The labeled components will then be allocated to
+% the individual runs, such that ICA exclusion can be performed on the run
+% data.
+
+% required functions:
+% - cm_label_ica_gui_20180116
+
+restoredefaultpath;
+clear all; close all; pack; clc;
+
+%% pathdef
+
+pn.eeg_root     = '/Volumes/LNDG/Projects/StateSwitch/dynamic/data/eeg/task/A_preproc/SA_preproc_study/';
+pn.EEG          = [pn.eeg_root, 'B_data/C_EEG_FT/'];
+pn.History      = [pn.eeg_root, 'B_data/D_History/'];
+% add ConMemEEG tools
+pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB/'];            addpath(genpath(pn.MWBtools));
+pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG/'];            addpath(genpath(pn.THGtools));
+pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common/'];         addpath(genpath(pn.commontools));
+pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+
+%% define Condition & IDs for preprocessing
+
+condEEG = 'dynamic';
+
+%% define IDs for visual screening
+
+ID = 1124;
+
+%% load data & start GUI
+
+% load data (full data)
+load([pn.EEG, num2str(ID), '_', condEEG, '_EEG_Rlm_Fhl_Ica.mat'],'data')
+
+if ~isfield(data.iclabels,'emg')
+    data.iclabels.emg = [];
+end
+
+% load config (ICA version)
+load([pn.History, num2str(ID), '_', condEEG, '_config.mat'],'config')
+
+% settings for ICA
+
+% ica config
+cfg.method           = 'runica';
+%cfg.channel          = {'all','-ECG','-A2'};
+cfg.trials           = 'all';
+cfg.numcomponent     = 'all';
+cfg.demean           = 'no';
+cfg.runica.extended  = 1;
+
+% ICA solution for segments
+cfg.unmixing     = data.unmixing;
+cfg.topolabel    = data.topolabel;
+
+% ICA (from weights)
+
+% components
+comp = ft_componentanalysis(cfg,data);
+    
+% include electrode information
+comp.elec = data.elec;
+
+% include ICA labeling
+comp.iclabels = data.iclabels;
+
+% clear cfg
+clear cfg
+
+% label ICs
+
+% settings
+cfg.topoall  = 'yes';
+cfg.chanlocs = data.chanlocs;
+
+% label electrodes
+data.iclabels = cm_label_ica_gui_20180116(cfg,comp);
+
+% keep labels in config
+config.ica1.iclabels.manual = data.iclabels;
+
+% clear variables
+clear cfg
+
+% save data
+
+% save data (full data)
+save([pn.EEG, num2str(ID), '_', condEEG, '_EEG_Rlm_Fhl_Ica.mat'],'data')
+
+% save config (ICA version)
+save([pn.History, num2str(ID), '_', condEEG, '_config.mat'],'config')
+
+% clear variables
+
+clear BLOCK ID comp config data

code/05_segmentation/G_STSW_segmentation_raw_data_180111.m

@@ -0,0 +1,224 @@
+function G_STSW_segmentation_raw_data_180111(id)
+
+%% G_STSW_segmentation_raw_data_180111
+
+% 170921 | JQK adapted from MD script
+% 180118 | adapted for STSW study YA, adapted for tardis
+
+%% initialize
+
+% restoredefaultpath;
+% clear all; close all; pack; clc;
+
+%% pathdef
+
+if ismac
+    pn.study        = '/Volumes/LNDG/Projects/StateSwitch/';
+    pn.dynamic_In   = [pn.eeg_root, 'B_data/B_EEG_ET_ByRun/'];
+    pn.triggerTiming= [pn.eeg_root, 'C_figures/D_TriggerTiming/'];
+    pn.EEG          = [pn.eeg_root, 'B_data/C_EEG_FT/']; mkdir(pn.EEG);
+    pn.History      = [pn.eeg_root, 'B_data/D_History/']; mkdir(pn.History);
+    % add ConMemEEG tools
+    pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB/'];           addpath(genpath(pn.MWBtools));
+    pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG/'];           addpath(genpath(pn.THGtools));
+    pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common/'];        addpath(genpath(pn.commontools));
+    pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+    pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+    pn.helper       = [pn.eeg_root, 'A_scripts/helper/'];           addpath(pn.helper);
+else
+    pn.root         = '/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/';
+    pn.EEG          = [pn.root, 'B_data/C_EEG_FT/'];
+    pn.History      = [pn.root, 'B_data/D_History/'];
+    pn.triggerTiming= [pn.root, 'C_figures/D_TriggerTiming/'];
+end
+
+%% define Condition & IDs for preprocessing
+
+condEEG = 'dynamic';
+
+%% define IDs for segmentation
+
+% N = 47;
+IDs = {'1117';'1118';'1120';'1124';'1126';'1131';'1132';'1135';'1136';'1138';'1144';'1151';'1158';'1160';'1163';'1164';'1167';'1169';'1172';'1173';'1178';'1182';'1215';'1216';'1219';'1221';'1223';'1227';'1228';'1233';'1234';'1237';'1239';'1240';'1243';'1245';'1247';'1250';'1252';'1257';'1261';'1265';'1266';'1268';'1270';'1276';'1281'};
+
+id = str2num(id);
+
+%  loop IDs
+clc;
+%for id = 1:length(IDs)
+    display(['processing ID ' num2str(IDs{id})]);
+    for iRun = 1:4
+        %%  load raw data
+
+        if ~exist([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_eyeEEG_Rlm_Fhl_rdSeg.mat'],'file')
+
+        % load config
+        load([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config.mat'],'config');
+
+        % copy ICA labeling info
+        configWithICA = load([pn.History, IDs{id}, '_', condEEG, '_config.mat'],'config');
+
+        config.ica1 = configWithICA.config.ica1;
+
+        %% ---- generate segmentation ---- %%
+
+        % load marker file
+        mrk = config.mrk;
+
+        for i = 1:size(mrk,2)
+            mrk_val{i,1} = mrk(1,i).value;
+        end
+
+        % generate trial structure
+        indOnset = find(strcmp(mrk_val(:,:),'S 17')); % (fix cue onset trigger = 'S 17')
+        indOnset = sortrows(indOnset,1);
+        indOffset = find(strcmp(mrk_val(:,:),'S 64')); % (ITI trigger = 'S 64')
+        indOffset = sortrows(indOffset,1);
+
+        h = figure;
+        plot(diff([mrk(indOnset).sample]))
+        hold on; plot(diff([mrk(indOffset).sample]))
+        pn.plotFolder = pn.triggerTiming;
+        figureName = ['A_triggerTiming_',IDs{id}, '_r',num2str(iRun)]; 
+        saveas(h, [pn.plotFolder, figureName], 'fig');
+        saveas(h, [pn.plotFolder, figureName], 'epsc');
+        saveas(h, [pn.plotFolder, figureName], 'png');
+        close(h);
+        
+        % Stim trials
+        trl = zeros(length(indOnset),3);
+        TOI1 = 1500; % segmentation before trigger
+        TOI2 = 1500; % segmentation following trigger
+
+        for j = 1:size(trl,1)
+            trl(j,1) = mrk(1,indOnset(j)).sample - TOI1; % segmentation from 1500 ms before stim
+            trl(j,2) = mrk(1,indOffset(j)).sample + TOI2; % to 1500 ms after trigger
+            trl(j,3) = -TOI1;% offset
+        end; clear j
+
+        % add trial structure to config
+        config.trl = trl;
+
+        % save config
+        save([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config.mat'],'config')
+
+        %% ----  load, filter, and re-reference raw data ---- %%
+
+        % define reading & preprocessing parameters
+        % read in all data first
+
+        cfg             = [];
+        %cfg.datafile    = [config.data_file];
+        cfg.datafile    = [pn.root, '/B_data/B_EEG_ET_ByRun/',IDs{id}, '_r',num2str(iRun), '_', condEEG, '_eyeEEG.set'];
+        cfg.trl         = config.trl;
+
+        cfg.channel     = {'all'};
+        cfg.implicitref = 'REF'; % recover implicit reference
+
+        data            = ft_preprocessing(cfg);
+
+         %% SWITCH CHANNELS ACCORDING TO ARRANGEMENT!
+
+        if max(strcmp(IDs{id}, {'1126'}))==1
+            data = SS_switchChannels_GreenYellow(data);       % green and yellow boxes exchanged
+            data = SS_switchChannels_Study_noA1(data);        % TP9 and TP10 exchanged manually
+        elseif max(strcmp(IDs{id}, {'1216'}))==1
+            data = SS_switchChannels_GreenYellow(data);       % green and yellow boxes exchanged
+            data = SS_switchChannels_Study(data);             % TP9(A1) and TP10(FCz) wrongly ordered in workspace
+        elseif max(strcmp(IDs{id}, {'1118'; '1215'; '1124'}))==1
+            data = SS_switchChannels_Study_noA1(data);        % TP9 and TP10 exchanged manually
+        else 
+            data = SS_switchChannels_Study(data);             % TP9(A1) and TP10(FCz) wrongly ordered in workspace
+        end
+        data = SS_changeChannels(data);
+
+        %%
+        % select eye data
+        cfg         =  [];
+        cfg.channel =  {'TIME', 'L_GAZE_X', 'L_GAZE_Y', ...
+            'L_AREA', 'L_VEL_X', 'L_VEL_Y', 'RES_X', 'RES_Y', 'INPUT'};
+
+        data_eye        = ft_preprocessing(cfg, data);
+        %%
+        % select EEG (non-eye/ECG/trigger) data
+        cfg                  =  [];
+        cfg.channel          = {'all', ...
+            '-TIME', '-L_GAZE_X', '-L_GAZE_Y', ...
+            '-L_AREA', '-L_VEL_X', '-L_VEL_Y', '-RES_X', '-RES_Y', '-INPUT', '-ECG'};
+
+        data_EEG             = ft_preprocessing(cfg, data);
+        %%
+        % filter & reref EEG data
+        cfg             = [];
+
+        cfg.continuous  = 'yes';
+        cfg.demean      = 'yes';
+
+        cfg.reref       = 'yes';
+        cfg.refchannel  = {'A1','A2'};
+        cfg.implicitref = 'A2';
+
+        cfg.hpfilter    = 'yes';
+        cfg.hpfreq      = .2;
+        cfg.hpfiltord   = 4;
+        cfg.hpfilttype  = 'but';
+
+        cfg.lpfilter    = 'yes';
+        cfg.lpfreq      = 125;
+        cfg.lpfiltord   = 4;
+        cfg.lpfilttype  = 'but';
+
+        % get data
+        data_EEG = ft_preprocessing(cfg, data_EEG);
+        flt = cfg;
+
+        % clear cfg structure
+        clear cfg
+
+        %% ----  resampling [500 Hz] ---- %%
+
+        % define settings for resampling
+        cfg.resamplefs = 500;                           % raw_eeg (ICA loop) at srate of 250Hz! now resample raw files (original files) to 500 Hz!
+        cfg.detrend    = 'no';
+        cfg.feedback   = 'no';
+        cfg.trials     = 'all';
+
+        % resample ALL data
+        data        = ft_resampledata(cfg,data);
+        if str2num(IDs{id})==1223 & (iRun == 2 | iRun == 3 | iRun == 4)
+            disp('Eye data skipped');
+        elseif str2num(IDs{id})==1228 & (iRun == 4)
+            disp('Eye data skipped');
+        else
+            data_eye    = ft_resampledata(cfg,data_eye);
+        end
+        data_EEG    = ft_resampledata(cfg,data_EEG);
+
+        resample = cfg;
+        % clear variables
+        clear cfg
+
+        % update config
+        config.seg.flt = flt;
+        config.seg.resample = resample;
+
+        % save config, data
+        save([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_eyeEEG_Rlm_Fhl_rdSeg'],'data')
+        if str2num(IDs{id})==1223 & (iRun == 2 | iRun == 3 | iRun == 4)
+            disp('Eye data skipped');
+        elseif str2num(IDs{id})==1228 & (iRun == 4)
+            disp('Eye data skipped');
+        else
+            save([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_eye_Rlm_Fhl_rdSeg'],'data_eye')
+        end
+        save([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_EEG_Rlm_Fhl_rdSeg'],'data_EEG')
+        save([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config'],'config')
+
+        % clear variables
+        clear config data data_eye data_EEG indOnset indOffset mrk mrk_val trl TOI1 TOI2
+
+        else
+
+        end % exist
+    end; clear iRun
+%end; clear id

code/05_segmentation/G_STSW_segmentation_raw_data_180111_START.sh

@@ -0,0 +1,22 @@
+#!/bin/bash
+
+# call the BOSC analysis by session and participant
+cd /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/G_STSW_segmentation/
+
+IDs=(1117 1118 1120 1124 1125 1126 1131 1132 1135 1136 1151 1160 1164 1167 1169 1172 1173 1178 1182 1214 1215 1216 1219 1223 1227 1228 1233 1234 1237 1239 1240 1243 1245 1247 1250 1252 1257 1261 1265 1266 1268 1270 1276 1281 2142 2253 2254 2255)
+
+mkdir /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/G_STSW_segment
+
+for i in $(seq 1 ${#IDs[@]}); do
+	echo "#PBS -N STSW_eeg_segment_${i}" 									> job
+	echo "#PBS -l walltime=03:00:00" 										>> job
+	echo "#PBS -l mem=4gb" 													>> job
+	echo "#PBS -j oe" 														>> job
+	echo "#PBS -o /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/G_STSW_segment" >> job
+	echo "#PBS -m n" 														>> job
+	echo "#PBS -d ." 														>> job
+	echo "./G_STSW_segmentation_raw_data_180111_run.sh /opt/matlab/R2016b $i " 	>> job
+	qsub job
+	rm job
+done
+done

code/05_segmentation/G_STSW_segmentation_raw_data_180111_noTardis.m

@@ -0,0 +1,251 @@
+function G_STSW_segmentation_raw_data_180111_noTardis()
+
+%% G_STSW_segmentation_raw_data_180111
+
+% 170921 | JQK adapted from MD script
+% 180118 | adapted for STSW study YA, adapted for tardis
+
+% config.data_file still refers to the non-tardis location
+
+% 180124 | skip data without eye info from processing
+
+%% initialize
+
+% restoredefaultpath;
+% clear all; close all; pack; clc;
+
+%% pathdef
+
+if ismac
+    pn.study        = '/Volumes/LNDG/Projects/StateSwitch/';
+    pn.eeg_root     = [pn.study, 'dynamic/data/eeg/task/A_preproc/SA_preproc_study/'];
+    pn.dynamic_In   = [pn.eeg_root, 'B_data/B_EEG_ET_ByRun/'];
+    pn.triggerTiming= [pn.eeg_root, 'C_figures/D_TriggerTiming/'];
+    pn.EEG          = [pn.eeg_root, 'B_data/C_EEG_FT/']; mkdir(pn.EEG);
+    pn.History      = [pn.eeg_root, 'B_data/D_History/']; mkdir(pn.History);
+    % add ConMemEEG tools
+    pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB/'];           addpath(genpath(pn.MWBtools));
+    pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG/'];           addpath(genpath(pn.THGtools));
+    pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common/'];        addpath(genpath(pn.commontools));
+    pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+    pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+    pn.helper       = [pn.eeg_root, 'A_scripts/helper/'];           addpath(pn.helper);
+else
+    pn.root         = '/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study/';
+    pn.EEG          = [pn.root, 'B_data/C_EEG_FT/'];
+    pn.History      = [pn.root, 'B_data/D_History/'];
+    pn.triggerTiming= [pn.root, 'C_figures/D_TriggerTiming/'];
+end
+
+%% define Condition & IDs for preprocessing
+
+condEEG = 'dynamic';
+
+%% define IDs for segmentation
+
+% N = 47 YAs + 53 OAs;
+IDs = {'1117';'1118';'1120';'1124';'1126';'1131';'1132';'1135';'1136';'1138';...
+    '1144';'1151';'1158';'1160';'1163';'1164';'1167';'1169';'1172';'1173';...
+    '1178';'1182';'1215';'1216';'1219';'1221';'1223';'1227';'1228';'1233';...
+    '1234';'1237';'1239';'1240';'1243';'1245';'1247';'1250';'1252';'1257';...
+    '1261';'1265';'1266';'1268';'1270';'1276';'1281';...
+    '2104';'2107';'2108';'2112';'2118';'2120';'2121';'2123';'2125';'2129';...
+    '2130';'2131';'2132';'2133';'2134';'2135';'2139';'2140';'2145';'2147';...
+    '2149';'2157';'2160';'2201';'2202';'2203';'2205';'2206';'2209';'2210';...
+    '2211';'2213';'2214';'2215';'2216';'2217';'2219';'2222';'2224';'2226';...
+    '2227';'2236';'2237';'2238';'2241';'2244';'2246';'2248';'2250';'2251';...
+    '2252';'2258';'2261'};
+
+%%  loop IDs
+
+clc;
+for id = 1:length(IDs)
+    display(['processing ID ' num2str(IDs{id})]);
+    for iRun = 1:4
+        %%  load raw data
+
+        if 1%~exist([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_eyeEEG_Rlm_Fhl_rdSeg.mat'],'file')
+
+        % load config
+        load([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config.mat'],'config');
+
+        % copy ICA labeling info
+        configWithICA = load([pn.History, IDs{id}, '_', condEEG, '_config.mat'],'config');
+
+        config.ica1 = configWithICA.config.ica1;
+
+        %% ---- generate segmentation ---- %%
+
+        % load marker file
+        mrk = config.mrk;
+
+        for i = 1:size(mrk,2)
+            mrk_val{i,1} = mrk(1,i).value;
+        end
+
+        % for 2160, the onset marker was not recorded: extract at beginning
+        if strcmp(IDs{id}, '2160') && iRun ==1
+            mrk_val{2,1} = 'S 17';
+        end
+        
+        % for 2203, the onset marker was not recorded: extract at beginning
+        if strcmp(IDs{id}, '2203') && iRun ==1
+            mrk_val{2,1} = 'S 17';
+        end
+        
+        % generate trial structure
+        indOnset = find(strcmp(mrk_val(:,:),'S 17')); % (fix cue onset trigger = 'S 17')
+        indOnset = sortrows(indOnset,1);
+        indOffset = find(strcmp(mrk_val(:,:),'S 64')); % (ITI trigger = 'S 64')
+        indOffset = sortrows(indOffset,1);
+
+        h = figure;
+        plot(diff([mrk(indOnset).sample]))
+        hold on; plot(diff([mrk(indOffset).sample]))
+        pn.plotFolder = pn.triggerTiming;
+        figureName = ['A_triggerTiming_',IDs{id}, '_r',num2str(iRun)]; 
+        saveas(h, [pn.plotFolder, figureName], 'png');
+        close(h);
+                
+        % Stim trials
+        trl = zeros(length(indOnset),3);
+        TOI1 = 1500; % segmentation before trigger
+        TOI2 = 1500; % segmentation following trigger
+
+        for j = 1:size(trl,1)
+            trl(j,1) = mrk(1,indOnset(j)).sample - TOI1; % segmentation from 1500 ms before stim
+            trl(j,2) = mrk(1,indOffset(j)).sample + TOI2; % to 1500 ms after trigger
+            trl(j,3) = -TOI1;% offset
+        end; clear j
+        
+        % if an onset occurs before recording onset: set to extract at 1
+        if numel(find(trl(:,1)<0))>0
+            trl(find(trl(:,1)<0),3) = trl(find(trl(:,1)<0),3)-trl(find(trl(:,1)<0),1);
+            trl(find(trl(:,1)<0),1) = 1;
+        end
+
+        % add trial structure to config
+        config.trl = trl;
+
+        % save config
+        save([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config.mat'],'config')
+
+        %% ----  load, filter, and re-reference raw data ---- %%
+
+        % define reading & preprocessing parameters
+        % read in all data first
+
+        cfg             = [];
+        cfg.datafile    = [config.data_file];                
+        cfg.trl         = config.trl;
+
+        cfg.channel     = {'all'};
+        cfg.implicitref = 'REF'; % recover implicit reference
+
+        data            = ft_preprocessing(cfg);
+
+         %% SWITCH CHANNELS ACCORDING TO ARRANGEMENT!
+
+        if max(strcmp(IDs{id}, {'1126'; '2227'}))==1
+            data = SS_switchChannels_GreenYellow(data);       % green and yellow boxes exchanged
+            data = SS_switchChannels_Study_noA1(data);        % TP9 and TP10 exchanged manually
+        elseif max(strcmp(IDs{id}, {'1216'}))==1
+            data = SS_switchChannels_GreenYellow(data);       % green and yellow boxes exchanged
+            data = SS_switchChannels_Study(data);             % TP9(A1) and TP10(FCz) wrongly ordered in workspace
+        elseif max(strcmp(IDs{id}, {'1118'; '1215'; '1124'}))==1
+            data = SS_switchChannels_Study_noA1(data);        % TP9 and TP10 exchanged manually
+        else 
+            data = SS_switchChannels_Study(data);             % TP9(A1) and TP10(FCz) wrongly ordered in workspace
+        end
+        data = SS_changeChannels(data);
+
+        %%
+        % select eye data
+        cfg         =  [];
+        cfg.channel =  {'TIME', 'L_GAZE_X', 'L_GAZE_Y', ...
+            'L_AREA', 'L_VEL_X', 'L_VEL_Y', 'RES_X', 'RES_Y', 'INPUT'};
+
+        data_eye        = ft_preprocessing(cfg, data);
+        %%
+        % select EEG (non-eye/ECG/trigger) data
+        cfg                  =  [];
+        cfg.channel          = {'all', ...
+            '-TIME', '-L_GAZE_X', '-L_GAZE_Y', ...
+            '-L_AREA', '-L_VEL_X', '-L_VEL_Y', '-RES_X', '-RES_Y', '-INPUT', '-ECG'};
+
+        data_EEG             = ft_preprocessing(cfg, data);
+        %%
+        % filter & reref EEG data
+        cfg             = [];
+
+        cfg.continuous  = 'yes';
+        cfg.demean      = 'yes';
+
+        cfg.reref       = 'yes';
+        cfg.refchannel  = {'A1','A2'};
+        cfg.implicitref = 'A2';
+
+        cfg.hpfilter    = 'yes';
+        cfg.hpfreq      = .2;
+        cfg.hpfiltord   = 4;
+        cfg.hpfilttype  = 'but';
+
+        cfg.lpfilter    = 'yes';
+        cfg.lpfreq      = 125;
+        cfg.lpfiltord   = 4;
+        cfg.lpfilttype  = 'but';
+
+        % get data
+        data_EEG = ft_preprocessing(cfg, data_EEG);
+        flt = cfg;
+
+        % clear cfg structure
+        clear cfg
+
+        %% ----  resampling [500 Hz] ---- %%
+
+        % define settings for resampling
+        cfg.resamplefs = 500;                           % raw_eeg (ICA loop) at srate of 250Hz! now resample raw files (original files) to 500 Hz!
+        cfg.detrend    = 'no';
+        cfg.feedback   = 'no';
+        cfg.trials     = 'all';
+
+        % resample ALL data
+        data        = ft_resampledata(cfg,data);
+        if str2num(IDs{id})==1223 & (iRun == 2 | iRun == 3 | iRun == 4)
+            disp('Eye data skipped');
+        elseif str2num(IDs{id})==1228 & (iRun == 4)
+            disp('Eye data skipped');
+        else
+            data_eye = ft_resampledata(cfg,data_eye);
+        end
+        data_EEG = ft_resampledata(cfg,data_EEG);
+
+        resample = cfg;
+        % clear variables
+        clear cfg
+
+        % update config
+        config.seg.flt = flt;
+        config.seg.resample = resample;
+
+        % save config, data
+        save([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_eyeEEG_Rlm_Fhl_rdSeg'],'data')
+        if str2num(IDs{id})==1223 & (iRun == 2 | iRun == 3 | iRun == 4)
+            disp('Eye data skipped');
+        elseif str2num(IDs{id})==1228 & (iRun == 4)
+            disp('Eye data skipped');
+        else
+            save([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_eye_Rlm_Fhl_rdSeg'],'data_eye')
+        end
+        save([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_EEG_Rlm_Fhl_rdSeg'],'data_EEG')
+        save([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config'],'config')
+
+        % clear variables
+        clear config data data_eye data_EEG indOnset indOffset mrk mrk_val trl TOI1 TOI2
+
+        else
+            disp(['Data already exists: ', IDs{id}, ' Run ',num2str(iRun)]);
+        end % exist
+    end; clear iRun
+end; clear id

code/05_segmentation/G_STSW_segmentation_raw_data_180111_prepare.sh

@@ -0,0 +1,21 @@
+#!/bin/bash
+
+# This script prepares tardis by compiling the necessary function in MATLAB.
+
+#ssh tardis # access tardis
+
+# check and choose matlab version
+#module avail matlab
+module load matlab/R2016b
+
+# compile functions
+
+matlab
+%% add fieldtrip toolbox
+addpath('/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/T_tools/fieldtrip-20170904/')
+ft_defaults()
+ft_compile_mex(true)
+%% go to analysis directory containing .m-file
+cd('/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/G_STSW_segmentation/')
+%% compile function and append dependencies
+mcc -m G_STSW_segmentation_raw_data_180111.m -a /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/helper

code/05_segmentation/G_STSW_segmentation_raw_data_180111_run.sh

@@ -0,0 +1,33 @@
+#!/bin/sh
+# script for execution of deployed applications
+#
+# Sets up the MATLAB Runtime environment for the current $ARCH and executes 
+# the specified command.
+#
+exe_name=$0
+exe_dir=`dirname "$0"`
+echo "------------------------------------------"
+if [ "x$1" = "x" ]; then
+  echo Usage:
+  echo    $0 \<deployedMCRroot\> args
+else
+  echo Setting up environment variables
+  MCRROOT="$1"
+  echo ---
+  LD_LIBRARY_PATH=.:${MCRROOT}/runtime/glnxa64 ;
+  LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/bin/glnxa64 ;
+  LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/sys/os/glnxa64;
+  LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/sys/opengl/lib/glnxa64;
+  export LD_LIBRARY_PATH;
+  echo LD_LIBRARY_PATH is ${LD_LIBRARY_PATH};
+  shift 1
+  args=
+  while [ $# -gt 0 ]; do
+      token=$1
+      args="${args} \"${token}\"" 
+      shift
+  done
+  eval "\"${exe_dir}/G_STSW_segmentation_raw_data_180111\"" $args
+fi
+exit
+

code/05_segmentation/mccExcludedFiles.log

@@ -0,0 +1,7 @@
+The List of Excluded Files
+Excluded files	 Exclusion Message ID	 Reason For Exclusion	 Exclusion Rule
+/opt/matlab/R2016b/toolbox/local/+matlab/+internal/+toolboxes/addInstalledToolboxesToJavaPath.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/[+]matlab/[+]internal/[+]toolboxes/
+/opt/matlab/R2016b/toolbox/local/+matlab/+internal/+toolboxes/addInstalledToolboxesToPath.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/[+]matlab/[+]internal/[+]toolboxes/
+/opt/matlab/R2016b/toolbox/local/+matlab/+internal/+zipAddOns/addInstalledZipAddOnsToPath.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/[+]matlab/[+]internal/[+]zipAddOns/
+/opt/matlab/R2016b/toolbox/local/pathdef.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/pathdef[.]m
+/opt/matlab/R2016b/toolbox/matlab/codetools/initdesktoputils.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/matlab/codetools/

code/06_automatic_artifact_correction/H_STSW_automatic_artifact_correction_180123.m

@@ -0,0 +1,278 @@
+function H_STSW_automatic_artifact_correction_180123(id)
+
+%% H_SS_automatic_artifact_correction_20170922
+
+%% initialize
+
+%restoredefaultpath;
+%clear all; close all; pack; clc;
+
+%% pathdef
+
+if ismac
+%     pn.study        = '/Volumes/LNDG/Projects/StateSwitch/';
+%     pn.eeg_root     = [pn.study, 'dynamic/data/eeg/task/A_preproc/SA_preproc_study_YA/'];
+    pn.eeg_root     = '/Users/kosciessa/Desktop/mountpoint_tardis_LNDG/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/';
+    pn.dynamic_In   = [pn.eeg_root, 'B_data/B_EEG_ET_ByRun/'];
+    pn.triggerTiming= [pn.eeg_root, 'C_figures/D_TriggerTiming/'];
+    pn.EEG          = [pn.eeg_root, 'B_data/C_EEG_FT/']; mkdir(pn.EEG);
+    pn.History      = [pn.eeg_root, 'B_data/D_History/']; mkdir(pn.History);
+    pn.logs         = [pn.eeg_root, 'Y_logs/H_ArtCorr/'];
+    % add ConMemEEG tools
+    pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB/'];           addpath(genpath(pn.MWBtools));
+    pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG/'];           addpath(genpath(pn.THGtools));
+    pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common/'];        addpath(genpath(pn.commontools));
+    pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+    pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+    pn.helper       = [pn.eeg_root, 'A_scripts/helper/'];           addpath(pn.helper);
+else
+    pn.root         = '/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/';
+    pn.EEG          = [pn.root, 'B_data/C_EEG_FT/'];
+    pn.History      = [pn.root, 'B_data/D_History/'];
+    pn.triggerTiming= [pn.root, 'C_figures/D_TriggerTiming/'];
+    pn.logs         = [pn.root, 'Y_logs/H_ArtCorr/'];
+end
+
+%% define Condition & IDs for preprocessing
+
+condEEG = 'dynamic';
+
+%% define IDs for segmentation
+
+% N = 47;
+IDs = {'1117';'1118';'1120';'1124';'1126';'1131';'1132';'1135';'1136';'1138';'1144';'1151';'1158';'1160';'1163';'1164';'1167';'1169';'1172';'1173';'1178';'1182';'1215';'1216';'1219';'1221';'1223';'1227';'1228';'1233';'1234';'1237';'1239';'1240';'1243';'1245';'1247';'1250';'1252';'1257';'1261';'1265';'1266';'1268';'1270';'1276';'1281'};
+
+id = str2num(id);
+
+%%  loop IDs
+
+ID_unavailable = cell(length(IDs),1);
+%for id = 1:length(IDs)
+    try
+        display(['processing ID ' num2str(IDs{id})]);
+        for iRun = 1:4
+        
+            % load data
+            load([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_EEG_Rlm_Fhl_rdSeg'],'data_EEG');
+            
+            dataICA = load([pn.EEG, IDs{id}, '_', condEEG, '_EEG_Rlm_Fhl_Ica'],'data');
+
+            data = data_EEG; clear data_EEG;
+            data.elec = dataICA.data.elec;
+            data.chanlocs = dataICA.data.chanlocs; clear dataICA;
+            
+           %% ------------------ ARTIFACT DETECTION - PREPARATION ----------------- %%
+
+            % load config
+
+            load([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config'],'config')
+            
+            config.elec = data.elec;
+
+            %% ICA (from weights)
+
+            % ica config
+            cfg.method           = 'runica';
+            cfg.channel          = {'all','-ECG','-A2'};
+            cfg.trials           = 'all';
+            cfg.numcomponent     = 'all';
+            cfg.demean           = 'yes';
+            cfg.runica.extended  = 1;
+
+            % ICA solution
+            cfg.unmixing     = config.ica1.unmixing;
+            cfg.topolabel    = config.ica1.topolabel;
+
+            % components
+            comp = ft_componentanalysis(cfg,data);
+
+            % clear cfg
+            clear cfg data
+
+            %% remove components
+
+            % get IC labels
+            iclabels = config.ica1.iclabels.manual;
+
+            % cfg for rejecting components (reject: blinks, eye movements,
+            % ecg, ref); JQK: added artifacts and emg
+            cfg.component = sortrows([iclabels.bli(:); iclabels.mov(:); iclabels.hrt(:); iclabels.ref(:); iclabels.art(:); iclabels.emg(:)])';
+            cfg.demean    = 'yes';
+
+            % reject components
+            data = ft_rejectcomponent(cfg,comp);
+
+            % clear cfg
+            clear cfg comp
+
+            %% remove eye & reference channels
+
+            cfg.channel     = {'all','-IOR','-LHEOG','-RHEOG','-A1'};
+            cfg.demean      = 'yes';
+
+            % remove channels
+            tmpdat = ft_preprocessing(cfg,data);
+
+            % clear cfg & data variable
+            clear cfg data
+
+
+            %% ------------------------- ARTIFACT DETECTION ------------------------ %%
+
+            % open log file
+            fid = fopen([pn.logs, 'log_' IDs{id}, '_r',num2str(iRun), '_' condEEG '_ArtCorr.txt'],'a');
+
+            % write log
+          % fprintf(fid,['*********************  ' BLOCK '  *********************\n']);                         % Undefined function or variable 'BLOCK'. but see: %% H_prep_data_for_analysis_20141218
+
+            fprintf(fid,['*********************   BLOCK   *********************\n']);
+            fprintf(fid,['function: H_SS_automatic_artifact_correction_20170922.m \n\n']);                
+            fprintf(fid,['eeg file = ' config.data_file '\n\n']);
+
+            n_trl = length(tmpdat.trial);
+
+            %%  get artifact contaminated channels by kurtosis, low & high frequency artifacts
+
+            cfg.criterion = 3;
+            cfg.recursive = 'no';
+
+            [index0 parm0 zval0] = cm_MWB_channel_x_epoch_artifacts_20170922(cfg,tmpdat);
+
+            % write log
+            tmp_log = '';
+            for j = 1:length(index0.c)
+                tmp_log = [tmp_log num2str(index0.c(j)) ' '];
+            end; clear j
+            tmp_log = [tmp_log(1:end-1) '\n'];
+            fprintf(fid,'(1) automatic bad channel detection:\n');
+            fprintf(fid,['MWB:          channel(s) ' tmp_log]);
+
+            % clear cfg
+            clear cfg tmp_log
+
+            %%  get artifact contaminated channels by FASTER
+
+            cfg.criterion = 3;
+            cfg.recursive = 'no';
+
+            [index1 parm1 zval1] = THG_FASTER_1_channel_artifacts_20140302(cfg,tmpdat);
+
+            % write log
+            tmp_log = '';
+            for j = 1:length(index1)
+                tmp_log = [tmp_log num2str(index1(j)) ' '];
+            end; clear j
+            tmp_log = [tmp_log(1:end-1) '\n'];
+            fprintf(fid,['FASTER:      channel(s) ' tmp_log]);
+
+            % clear cfg
+            clear cfg tmp_log
+
+            %%  interpolate artifact contaminated channels
+
+            % collect bad channels
+            badchan = unique([index0.c; index1]);
+
+            fprintf(fid,['--> ' num2str(length(badchan)) ' channels interpolated\n\n']);
+            
+            cfg.method     = 'spline';
+            cfg.badchannel = tmpdat.label(badchan);
+            cfg.trials     = 'all';
+            cfg.lambda     = 1e-5; 
+            cfg.order      = 4; 
+            cfg.elec       = config.elec;
+
+            % interpolate
+            tmpdat = ft_channelrepair(cfg,tmpdat);
+
+            % clear cfg
+            clear cfg
+
+            %% equalize duration of trials to the trial with fewest samples
+            
+            for indTrial = 1:64
+                TrialLength(indTrial) = size(tmpdat.trial{indTrial},2);
+            end
+            [val, idx] = min(TrialLength);
+            for indTrial = 1:64
+                tmpdat.trial{indTrial} = tmpdat.trial{indTrial}(:,1:val);
+                tmpdat.time{indTrial} = tmpdat.time{idx};
+            end
+            
+            %%  get artifact contaminated epochs & exclude epochs
+            % includes: - THG_MWB_channel_x_epoch_artifacts_20140311
+            %           - THG_FASTER_2_epoch_artifacts_20140302
+            % recursive epoch exclusion!
+
+            [tmpdat index] = THG_automatic_artifact_correction_trials_20170922(tmpdat);
+
+            % write log
+            fprintf(fid,'(2) automatic recursive bad epoch detection:\n');
+            fprintf(fid,['MWB & FASTER: ' num2str(n_trl-length(index)) ' bad epoch(s) detected\n\n']);
+
+            %%  get channel x epoch artifacts
+
+            % cfg
+            cfg.criterion = 3;
+            cfg.recursive = 'yes';
+
+            [index3 parm3 zval3] = THG_FASTER_4_channel_x_epoch_artifacts_20140302(cfg,tmpdat);
+
+            % write log
+            fprintf(fid,'(3) automatic single epoch/channel detection:\n');
+            fprintf(fid,['FASTER:       ' num2str(sum(sum(index3))) ' channel(s) x trial(s) detected\n\n']);
+
+            % clear cfg
+            clear cfg
+
+            %%  collect and save detected artifacts & artifact correction infos
+
+            % include ArtDect.parameters
+
+            % bad channels
+            ArtDect.channels = badchan;
+
+            % bad trials
+            tmp  = zeros(n_trl,1); tmp(index,1) = 1;
+            badtrl = find(tmp==0);
+            ArtDect.trials  = badtrl; 
+            clear tmp
+
+            % bad single epoch(s)/channel(s) - after exclusion of bad epochs
+            ArtDect.channels_x_trials = index3;
+
+            % overview
+            ind = [1:n_trl];
+            ind(badtrl) = [];
+            tmp = ones(length(tmpdat.label),n_trl);
+            tmp(:,ind) = index3;
+            tmp(badchan,:) = 1;
+            ArtDect.channels_x_trials_all = tmp;
+            clear tmp
+
+            % save version
+            ArtDect.version = 'JQK 20180123';
+
+            % add to config
+            config.ArtDect = ArtDect;
+
+            % save config
+            save([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config'],'config')
+            
+            %%  finalize log
+
+            % log
+            fprintf(fid,'Artifact detection completed.\n');
+            fprintf(fid,['Information saved to: ' IDs{id} '_config.mat\n\n']);
+            fclose(fid);
+        end % run
+    catch ME
+        warning('Error occured. Please check.');
+        rethrow(ME)
+        fprintf('\n ID not availble! Skip! \n')
+        ID_unavailable{id,1} = (IDs{id});
+    end
+    
+%end; clear id
+
+end

code/06_automatic_artifact_correction/H_STSW_automatic_artifact_correction_180123_START.sh

@@ -0,0 +1,22 @@
+#!/bin/bash
+
+# call the BOSC analysis by session and participant
+cd /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/H_STSW_automatic_artifact_correction/
+
+IDs=(1117 1118 1120 1124 1125 1126 1131 1132 1135 1136 1151 1160 1164 1167 1169 1172 1173 1178 1182 1214 1215 1216 1219 1223 1227 1228 1233 1234 1237 1239 1240 1243 1245 1247 1250 1252 1257 1261 1265 1266 1268 1270 1276 1281 2142 2253 2254 2255)
+
+mkdir /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/H_ArtCorr
+
+for i in $(seq 1 ${#IDs[@]}); do
+	echo "#PBS -N STSW_eeg_artCor_${i}" 									> job
+	echo "#PBS -l walltime=04:00:00" 										>> job
+	echo "#PBS -l mem=8gb" 													>> job
+	echo "#PBS -j oe" 														>> job
+	echo "#PBS -o /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/H_ArtCorr" >> job
+	echo "#PBS -m n" 														>> job
+	echo "#PBS -d ." 														>> job
+	echo "./H_STSW_automatic_artifact_correction_180123_run.sh /opt/matlab/R2016b $i " 	>> job
+	qsub job
+	rm job
+done
+done

code/06_automatic_artifact_correction/H_STSW_automatic_artifact_correction_180123_START_rerun.sh

@@ -0,0 +1,30 @@
+#!/bin/bash
+
+# call the BOSC analysis by session and participant
+cd /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/H_STSW_automatic_artifact_correction/
+
+IDs=(1117 1118 1120 1124 1125 1126 1131 1132 1135 1136 1151 1160 1164 1167 1169 1172 1173 1178 1182 1214 1215 1216 1219 1223 1227 1228 1233 1234 1237 1239 1240 1243 1245 1247 1250 1252 1257 1261 1265 1266 1268 1270 1276 1281 2142 2253 2254 2255)
+
+mkdir /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/H_ArtCorr
+
+	echo "#PBS -N STSW_eeg_artCor_27" 									> job
+	echo "#PBS -l walltime=04:00:00" 										>> job
+	echo "#PBS -l mem=8gb" 													>> job
+	echo "#PBS -j oe" 														>> job
+	echo "#PBS -o /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/H_ArtCorr" >> job
+	echo "#PBS -m n" 														>> job
+	echo "#PBS -d ." 														>> job
+	echo "./H_STSW_automatic_artifact_correction_180123_run.sh /opt/matlab/R2016b 27 " 	>> job
+	qsub job
+	rm job
+
+	echo "#PBS -N STSW_eeg_artCor_29" 									> job
+	echo "#PBS -l walltime=04:00:00" 										>> job
+	echo "#PBS -l mem=8gb" 													>> job
+	echo "#PBS -j oe" 														>> job
+	echo "#PBS -o /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/Y_logs/H_ArtCorr" >> job
+	echo "#PBS -m n" 														>> job
+	echo "#PBS -d ." 														>> job
+	echo "./H_STSW_automatic_artifact_correction_180123_run.sh /opt/matlab/R2016b 29 " 	>> job
+	qsub job
+	rm job

code/06_automatic_artifact_correction/H_STSW_automatic_artifact_correction_180123_noTardis.m

@@ -0,0 +1,307 @@
+function H_STSW_automatic_artifact_correction_180123_noTardis()
+
+%% H_SS_automatic_artifact_correction_20170922
+
+% 180326 | instead of equalizing to smallest amount of samples, interpolate
+%           for the same timing across subjects. Note that this does not
+%           counteract any natural jitter in the paradigm and should only
+%           equalize the size of matrices (fill with zeros! at ends).
+%           Artefact detection does NOT currently work with potential NaN
+%           values that would naturally arise from the interpolation.
+
+%% initialize
+
+%restoredefaultpath;
+%clear all; close all; pack; clc;
+
+%% pathdef
+
+if ismac
+    pn.study        = '/Volumes/LNDG/Projects/StateSwitch/';
+    pn.eeg_root     = [pn.study, 'dynamic/data/eeg/task/A_preproc/SA_preproc_study/'];
+    pn.dynamic_In   = [pn.eeg_root, 'B_data/B_EEG_ET_ByRun/'];
+    pn.triggerTiming= [pn.eeg_root, 'C_figures/D_TriggerTiming/'];
+    pn.EEG          = [pn.eeg_root, 'B_data/C_EEG_FT/']; mkdir(pn.EEG);
+    pn.History      = [pn.eeg_root, 'B_data/D_History/']; mkdir(pn.History);
+    pn.logs         = [pn.eeg_root, 'Y_logs/H_ArtCorr/'];
+    % add ConMemEEG tools
+    pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB/'];           addpath(genpath(pn.MWBtools));
+    pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG/'];           addpath(genpath(pn.THGtools));
+    pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common/'];        addpath(genpath(pn.commontools));
+    pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+    pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+    pn.helper       = [pn.eeg_root, 'A_scripts/helper/'];           addpath(pn.helper);
+else
+    pn.root         = '/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study/';
+    pn.EEG          = [pn.root, 'B_data/C_EEG_FT/'];
+    pn.History      = [pn.root, 'B_data/D_History/'];
+    pn.triggerTiming= [pn.root, 'C_figures/D_TriggerTiming/'];
+    pn.logs         = [pn.root, 'Y_logs/H_ArtCorr/'];
+end
+
+%% define Condition & IDs for preprocessing
+
+condEEG = 'dynamic';
+
+%% define IDs for segmentation
+
+% N = 47 YAs + 53 OAs;
+IDs = {'1117';'1118';'1120';'1124';'1126';'1131';'1132';'1135';'1136';'1138';...
+    '1144';'1151';'1158';'1160';'1163';'1164';'1167';'1169';'1172';'1173';...
+    '1178';'1182';'1215';'1216';'1219';'1221';'1223';'1227';'1228';'1233';...
+    '1234';'1237';'1239';'1240';'1243';'1245';'1247';'1250';'1252';'1257';...
+    '1261';'1265';'1266';'1268';'1270';'1276';'1281';...
+    '2104';'2107';'2108';'2112';'2118';'2120';'2121';'2123';'2125';'2129';...
+    '2130';'2131';'2132';'2133';'2134';'2135';'2139';'2140';'2145';'2147';...
+    '2149';'2157';'2160';'2201';'2202';'2203';'2205';'2206';'2209';'2210';...
+    '2211';'2213';'2214';'2215';'2216';'2217';'2219';'2222';'2224';'2226';...
+    '2227';'2236';'2237';'2238';'2241';'2244';'2246';'2248';'2250';'2251';...
+    '2252';'2258';'2261'};
+
+%%  loop IDs
+
+ID_unavailable = cell(length(IDs),1);
+for id = 1:length(IDs)
+    try
+        display(['processing ID ' num2str(IDs{id})]);
+        for iRun = 1:4
+        
+            % load config
+            load([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config'],'config')
+            
+            if ~isfield(config, 'ArtDect')
+            
+                % load data
+                load([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_EEG_Rlm_Fhl_rdSeg'],'data_EEG');
+
+                dataICA = load([pn.EEG, IDs{id}, '_', condEEG, '_EEG_Rlm_Fhl_Ica'],'data');
+
+                data = data_EEG; clear data_EEG;
+                data.elec = dataICA.data.elec;
+                data.chanlocs = dataICA.data.chanlocs; clear dataICA;
+
+               %% ------------------ ARTIFACT DETECTION - PREPARATION ----------------- %%
+
+                config.elec = data.elec;
+
+                %% ICA (from weights)
+
+                % ica config
+                cfg.method           = 'runica';
+                cfg.channel          = {'all','-ECG','-A2'};
+                cfg.trials           = 'all';
+                cfg.numcomponent     = 'all';
+                cfg.demean           = 'yes';
+                cfg.runica.extended  = 1;
+
+                % ICA solution
+                cfg.unmixing     = config.ica1.unmixing;
+                cfg.topolabel    = config.ica1.topolabel;
+
+                % components
+                comp = ft_componentanalysis(cfg,data);
+
+                % clear cfg
+                clear cfg data
+
+                %% remove components
+
+                % get IC labels
+                iclabels = config.ica1.iclabels.manual;
+
+                % cfg for rejecting components (reject: blinks, eye movements,
+                % ecg, ref); JQK: added artifacts and emg
+                cfg.component = sortrows([iclabels.bli(:); iclabels.mov(:); iclabels.hrt(:); iclabels.ref(:); iclabels.art(:); iclabels.emg(:)])';
+                cfg.demean    = 'yes';
+
+                % reject components
+                data = ft_rejectcomponent(cfg,comp);
+
+                % clear cfg
+                clear cfg comp
+
+                %% remove eye & reference channels
+
+                cfg.channel     = {'all','-IOR','-LHEOG','-RHEOG','-A1'};
+                cfg.demean      = 'yes';
+
+                % remove channels
+                tmpdat = ft_preprocessing(cfg,data);
+
+                % clear cfg & data variable
+                clear cfg data
+
+
+                %% ------------------------- ARTIFACT DETECTION ------------------------ %%
+
+                % open log file
+                fid = fopen([pn.logs, 'log_' IDs{id}, '_r',num2str(iRun), '_' condEEG '_ArtCorr.txt'],'a');
+
+                % write log
+              % fprintf(fid,['*********************  ' BLOCK '  *********************\n']);                         % Undefined function or variable 'BLOCK'. but see: %% H_prep_data_for_analysis_20141218
+
+                fprintf(fid,['*********************   BLOCK   *********************\n']);
+                fprintf(fid,['function: H_SS_automatic_artifact_correction_20170922.m \n\n']);                
+                fprintf(fid,['eeg file = ' config.data_file '\n\n']);
+
+                n_trl = length(tmpdat.trial);
+
+                %%  get artifact contaminated channels by kurtosis, low & high frequency artifacts
+
+                cfg.criterion = 3;
+                cfg.recursive = 'no';
+
+                [index0 parm0 zval0] = cm_MWB_channel_x_epoch_artifacts_20170922(cfg,tmpdat);
+
+                % write log
+                tmp_log = '';
+                for j = 1:length(index0.c)
+                    tmp_log = [tmp_log num2str(index0.c(j)) ' '];
+                end; clear j
+                tmp_log = [tmp_log(1:end-1) '\n'];
+                fprintf(fid,'(1) automatic bad channel detection:\n');
+                fprintf(fid,['MWB:          channel(s) ' tmp_log]);
+
+                % clear cfg
+                clear cfg tmp_log
+
+                %%  get artifact contaminated channels by FASTER
+
+                cfg.criterion = 3;
+                cfg.recursive = 'no';
+
+                [index1 parm1 zval1] = THG_FASTER_1_channel_artifacts_20140302(cfg,tmpdat);
+
+                % write log
+                tmp_log = '';
+                for j = 1:length(index1)
+                    tmp_log = [tmp_log num2str(index1(j)) ' '];
+                end; clear j
+                tmp_log = [tmp_log(1:end-1) '\n'];
+                fprintf(fid,['FASTER:      channel(s) ' tmp_log]);
+
+                % clear cfg
+                clear cfg tmp_log
+
+                %%  interpolate artifact contaminated channels
+
+                % collect bad channels
+                badchan = unique([index0.c; index1]);
+
+                fprintf(fid,['--> ' num2str(length(badchan)) ' channels interpolated\n\n']);
+
+                cfg.method     = 'spline';
+                cfg.badchannel = tmpdat.label(badchan);
+                cfg.trials     = 'all';
+                cfg.lambda     = 1e-5; 
+                cfg.order      = 4; 
+                cfg.elec       = config.elec;
+
+                % interpolate
+                tmpdat = ft_channelrepair(cfg,tmpdat);
+
+                % clear cfg
+                clear cfg
+
+                %% interpolate data to requested time (necessary due to acquisition jitter)
+                
+                % 2160 has only 62 trials, i.e. the first two trials were
+                % not recorded.
+                
+                timeRequested = -1.5:1/500:9.566;
+                
+                for indTrial = 1:size(tmpdat.trial,2)
+                    recTime = tmpdat.time{indTrial};
+                    recData = tmpdat.trial{indTrial};
+                    for indChannel = 1:size(recData,1)
+                        interpData(indChannel,:) = interp1(recTime,recData(indChannel,:),timeRequested,'linear');
+                    end
+                    interpData(isnan(interpData)) = 0;
+                    tmpdat.trial{indTrial} = interpData;
+                    tmpdat.time{indTrial} = timeRequested;
+                end
+                
+                %figure; plot(timeRequested,interpData, 'k'); hold on; plot(recTime,recData(1,:), 'r')
+                
+                %%  get artifact contaminated epochs & exclude epochs
+                % includes: - THG_MWB_channel_x_epoch_artifacts_20140311
+                %           - THG_FASTER_2_epoch_artifacts_20140302
+                % recursive epoch exclusion!
+
+                [tmpdat index] = THG_automatic_artifact_correction_trials_20170922(tmpdat);
+
+                % write log
+                fprintf(fid,'(2) automatic recursive bad epoch detection:\n');
+                fprintf(fid,['MWB & FASTER: ' num2str(n_trl-length(index)) ' bad epoch(s) detected\n\n']);
+
+                %%  get channel x epoch artifacts
+
+                % cfg
+                cfg.criterion = 3;
+                cfg.recursive = 'yes';
+
+                [index3 parm3 zval3] = THG_FASTER_4_channel_x_epoch_artifacts_20140302(cfg,tmpdat);
+
+                % write log
+                fprintf(fid,'(3) automatic single epoch/channel detection:\n');
+                fprintf(fid,['FASTER:       ' num2str(sum(sum(index3))) ' channel(s) x trial(s) detected\n\n']);
+
+                % clear cfg
+                clear cfg
+
+                %%  collect and save detected artifacts & artifact correction infos
+
+                % include ArtDect.parameters
+
+                % bad channels
+                ArtDect.channels = badchan;
+
+                % bad trials
+                tmp  = zeros(n_trl,1); tmp(index,1) = 1;
+                badtrl = find(tmp==0);
+                ArtDect.trials  = badtrl; 
+                clear tmp
+
+                % bad single epoch(s)/channel(s) - after exclusion of bad epochs
+                ArtDect.channels_x_trials = index3;
+
+                % overview
+                ind = [1:n_trl];
+                ind(badtrl) = [];
+                tmp = ones(length(tmpdat.label),n_trl);
+                tmp(:,ind) = index3;
+                tmp(badchan,:) = 1;
+                ArtDect.channels_x_trials_all = tmp;
+                clear tmp
+
+                % save version
+                ArtDect.version = 'JQK 20180123';
+
+                % add to config
+                config.ArtDect = ArtDect;
+
+                % save config
+                save([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config'],'config')
+
+                %%  finalize log
+
+                % log
+                fprintf(fid,'Artifact detection completed.\n');
+                fprintf(fid,['Information saved to: ' IDs{id} '_config.mat\n\n']);
+                fclose(fid);
+                
+                
+            else
+                disp([IDs{id}, ' Run ',num2str(iRun), ' has already run through automatic artifact detection']);
+            end % if data was not processed yet
+        end % run
+    catch ME
+        warning('Error occured. Please check.');
+        rethrow(ME)
+        fprintf('\n ID not availble! Skip! \n')
+        ID_unavailable{id,1} = (IDs{id});
+    end
+    
+end; clear id
+
+end

code/06_automatic_artifact_correction/H_STSW_automatic_artifact_correction_180123_prepare.sh

@@ -0,0 +1,21 @@
+#!/bin/bash
+
+# This script prepares tardis by compiling the necessary function in MATLAB.
+
+#ssh tardis # access tardis
+
+# check and choose matlab version
+#module avail matlab
+module load matlab/R2016b
+
+# compile functions
+
+matlab
+%% add fieldtrip toolbox
+addpath('/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/T_tools/fieldtrip-20170904/')
+ft_defaults()
+ft_compile_mex(true)
+%% go to analysis directory containing .m-file
+cd('/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/A_scripts/H_STSW_automatic_artifact_correction/')
+%% compile function and append dependencies
+mcc -m H_STSW_automatic_artifact_correction_180123.m -a /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/T_tools/fnct_common -a /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/T_tools/fnct_THG -a /home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study_YA/T_tools/fnct_MWB/MWB_ArtfDetec

code/06_automatic_artifact_correction/H_STSW_automatic_artifact_correction_180123_run.sh

@@ -0,0 +1,33 @@
+#!/bin/sh
+# script for execution of deployed applications
+#
+# Sets up the MATLAB Runtime environment for the current $ARCH and executes 
+# the specified command.
+#
+exe_name=$0
+exe_dir=`dirname "$0"`
+echo "------------------------------------------"
+if [ "x$1" = "x" ]; then
+  echo Usage:
+  echo    $0 \<deployedMCRroot\> args
+else
+  echo Setting up environment variables
+  MCRROOT="$1"
+  echo ---
+  LD_LIBRARY_PATH=.:${MCRROOT}/runtime/glnxa64 ;
+  LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/bin/glnxa64 ;
+  LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/sys/os/glnxa64;
+  LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/sys/opengl/lib/glnxa64;
+  export LD_LIBRARY_PATH;
+  echo LD_LIBRARY_PATH is ${LD_LIBRARY_PATH};
+  shift 1
+  args=
+  while [ $# -gt 0 ]; do
+      token=$1
+      args="${args} \"${token}\"" 
+      shift
+  done
+  eval "\"${exe_dir}/H_STSW_automatic_artifact_correction_180123\"" $args
+fi
+exit
+

code/06_automatic_artifact_correction/mccExcludedFiles.log

@@ -0,0 +1,7 @@
+The List of Excluded Files
+Excluded files	 Exclusion Message ID	 Reason For Exclusion	 Exclusion Rule
+/opt/matlab/R2016b/toolbox/local/+matlab/+internal/+toolboxes/addInstalledToolboxesToJavaPath.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/[+]matlab/[+]internal/[+]toolboxes/
+/opt/matlab/R2016b/toolbox/local/+matlab/+internal/+toolboxes/addInstalledToolboxesToPath.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/[+]matlab/[+]internal/[+]toolboxes/
+/opt/matlab/R2016b/toolbox/local/+matlab/+internal/+zipAddOns/addInstalledZipAddOnsToPath.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/[+]matlab/[+]internal/[+]zipAddOns/
+/opt/matlab/R2016b/toolbox/local/pathdef.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/local/pathdef[.]m
+/opt/matlab/R2016b/toolbox/matlab/codetools/initdesktoputils.m	MATLAB:depfun:req:ExcludedBy	Cannot be packaged for use in the target environment MCR. Removed from the parts list by license Compiler.	/opt/matlab/R2016b/toolbox/matlab/codetools/

code/07_STSW_prep_data_for_analysis.m

@@ -0,0 +1,306 @@
+%% I_SS_prep_data_for_analysis_20170922
+
+% 180124 | now also exclude ART and EMG artifacts
+% 180327 | include OAs, deal wih dropped trials of subjects with missing onsets
+%        | interpolate data to requested times vs. choosing shortest
+%        duration (zero interpolation)
+
+% This script fixes issues with missing onsets in the following way:
+
+% IMPORTANT: For 2160, the first two trials were not available,
+% hence we have to add +2 to the index to identify the correct
+% behavioral trials that belong to the EEG data. For 2203 the
+% first trial is missing. The third and second trial should
+% also be removed, regardless of whether they are indicated to
+% be dropped due to automatic correction.
+
+%% initialize
+
+restoredefaultpath;
+clear all; close all; pack; clc;
+
+%% pathdef
+
+if ismac
+    pn.study        = '/Volumes/LNDG/Projects/StateSwitch/';
+    pn.eeg_root     = [pn.study, 'dynamic/data/eeg/task/A_preproc/SA_preproc_study/'];
+    pn.EEG          = [pn.eeg_root, 'B_data/C_EEG_FT/']; %mkdir(pn.EEG);
+    pn.History      = [pn.eeg_root, 'B_data/D_History/']; %mkdir(pn.History);
+    % add ConMemEEG tools
+    pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB/'];           addpath(genpath(pn.MWBtools));
+    pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG/'];           addpath(genpath(pn.THGtools));
+    pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common/'];        addpath(genpath(pn.commontools));
+    pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+    pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+    pn.helper       = [pn.eeg_root, 'A_scripts/helper/'];           addpath(pn.helper);
+else
+    pn.root         = '/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study/';
+    pn.EEG          = [pn.root, 'B_data/C_EEG_FT/'];
+    pn.History      = [pn.root, 'B_data/D_History/'];
+end
+
+%% define Condition & IDs for preprocessing
+
+condEEG = 'dynamic';
+
+%% define IDs for segmentation
+
+% N = 47 YAs + 53 OAs;
+IDs = {'1117';'1118';'1120';'1124';'1126';'1131';'1132';'1135';'1136';'1138';...
+    '1144';'1151';'1158';'1160';'1163';'1164';'1167';'1169';'1172';'1173';...
+    '1178';'1182';'1215';'1216';'1219';'1221';'1223';'1227';'1228';'1233';...
+    '1234';'1237';'1239';'1240';'1243';'1245';'1247';'1250';'1252';'1257';...
+    '1261';'1265';'1266';'1268';'1270';'1276';'1281';...
+    '2104';'2107';'2108';'2112';'2118';'2120';'2121';'2123';'2125';'2129';...
+    '2130';'2131';'2132';'2133';'2134';'2135';'2139';'2140';'2145';'2147';...
+    '2149';'2157';'2160';'2201';'2202';'2203';'2205';'2206';'2209';'2210';...
+    '2211';'2213';'2214';'2215';'2216';'2217';'2219';'2222';'2224';'2226';...
+    '2227';'2236';'2237';'2238';'2241';'2244';'2246';'2248';'2250';'2251';...
+    '2252';'2258';'2261'};
+
+%%  loop IDs
+
+setting.plot = 0; % plot imagesc of resulting timeseries
+
+% prior to loop create matrixes for overviews over artifact-free trials and total trial numbers
+overview_nartfree_trials = {};
+overview_trial_numbers = {};
+ID_unavailable = cell(length(IDs),4);
+for id = 1:length(IDs)
+    for iRun = 1:4
+        
+%         if ~exist([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_EEG_Rlm_Fhl_rdSeg_Art.mat'])
+%             display(['processing ID ' IDs{id}, ' Run: ', num2str(iRun)]);
+%         else
+%             display(['ID ' IDs{id}, ' Run: ', num2str(iRun), ' has already been processed.']);
+%             continue;
+%         end
+        
+        try
+        %%  load files
+            
+            % load data
+            load([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_EEG_Rlm_Fhl_rdSeg'],'data_EEG');
+            dataICA = load([pn.EEG, IDs{id}, '_', condEEG, '_EEG_Rlm_Fhl_Ica'],'data');
+            data = data_EEG; clear data_EEG;
+            data.elec = dataICA.data.elec;
+            data.chanlocs = dataICA.data.chanlocs; clear dataICA;
+            % load config
+            load([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config'],'config')
+
+        %% cut data to requested time (necessary due to acquisition jitter)
+        % Note that this does not affect the data, but we simply cut some
+        % sample points from the ITI which is not of interest anyways.
+
+            trig.eventDurInS = [1.5, 2, 1, 3, 2, 1.5]; % total time to extract is 11s
+            trig.eventTriggers = {'S 17'; 'S 72'; 'S 20'; 'S 24'; 'S 64'};
+        
+            % 2160 has only 62 trials, i.e. the first two trials were
+            % not recorded.
+
+            cfg = [];
+            cfg.toilim = [-1.5 9.5];
+            data = ft_redefinetrial(cfg, data);
+            
+        %%  define parameter
+
+            % ICs to reject
+            iclabels    = config.ica1.iclabels.manual;
+            %ics2reject = [iclabels.bli(:); iclabels.mov(:); iclabels.hrt(:); iclabels.ref(:)];
+            ics2reject = [iclabels.bli(:); iclabels.mov(:); iclabels.hrt(:); iclabels.ref(:); iclabels.art(:); iclabels.emg(:)];
+            clear iclabels
+
+            % trials to remove
+            trls2remove = config.ArtDect.trials;
+
+            % channels to interpolate
+            chns2interp = config.ArtDect.channels;
+
+        %% for select IDs, interpolate channels prior to ICA
+        
+        if strcmp(IDs{id}, '2112')
+            cfg = [];
+            cfg.channel = {'all','-ECG','-A2'};
+            dataSelect = ft_preprocessing(cfg,data);
+            cfg = [];
+            cfg.method     = 'spline';
+            if strcmp(IDs{id}, '2112')
+                cfg.badchannel = {'Fz'; 'FCz'};
+            end
+            cfg.trials     = 'all';
+            cfg.lambda     = 1e-5; 
+            cfg.order      = 4; 
+            data = ft_channelrepair(cfg,data);
+            clear cfg;
+        end
+            
+        %% ICA (from weights)
+
+            % ica config
+            cfg.method           = 'runica';
+            cfg.channel          = {'all','-ECG','-A2'};
+            cfg.trials           = 'all';
+            cfg.numcomponent     = 'all';
+            cfg.demean           = 'yes';
+            cfg.runica.extended  = 1;
+
+            % ICA solution
+            cfg.unmixing     = config.ica1.unmixing;
+            cfg.topolabel    = config.ica1.topolabel;
+
+            % components
+            comp = ft_componentanalysis(cfg,data);
+
+            % clear cfg
+            clear cfg data
+
+        %% remove components
+
+            % cfg for rejecting components (reject: blinks, eye movements, ecg, ref)
+            cfg.component = sortrows(ics2reject)';
+            cfg.demean    = 'yes';
+
+            % reject components
+            data = ft_rejectcomponent(cfg,comp);
+
+            % clear cfg
+            clear cfg comp ics2reject
+
+        %%  remove trials
+
+            % IMPORTANT: For 2160, the first two trials were not available,
+            % hence we have to add +2 to the index to identify the correct
+            % behavioral trials that belong to the EEG data. For 2203 the
+            % first trial is missing. The third and second trial should
+            % also be removed, regardless of whether they are indicated to
+            % be dropped due to automatic correction.
+        
+            if strcmp(IDs{id}, '2160') && iRun==1
+                trls2remove=[trls2remove;1];
+            elseif strcmp(IDs{id}, '2203') && iRun==1
+                trls2remove=[trls2remove;1];
+            end
+            
+            % define trials to keep
+            trials               = 1:length(data.trial);
+            trials(trls2remove)  = [];
+
+            if strcmp(IDs{id}, '2160') && iRun==1
+                trials=trials+2;
+                data.time = cat(2,data.time{end},data.time{end},data.time);
+                data.trial = cat(2,data.trial{end},data.trial{end},data.trial);
+            elseif strcmp(IDs{id}, '2203') && iRun==1
+                trials=trials+1;
+                data.time = cat(2,data.time{end},data.time);
+                data.trial = cat(2,data.trial{end},data.trial);
+            end
+            
+            % config for deleting trials
+            cfg.trials   = trials;
+
+            % remove trials
+            data   = ft_preprocessing(cfg,data);
+            nartfree_trials = trials;
+            overview_nartfree_trials{id,iRun} = nartfree_trials;
+
+            config.nartfree_trials = nartfree_trials;
+            
+            % clear variables
+            clear cfg trials trls2remove
+
+        %% remove eye & reference channels (before interpolation)
+
+            cfg.channel     = {'all','-IOR','-LHEOG','-RHEOG','-A1'};
+            cfg.demean      = 'yes';
+
+            % remove channels
+            tmpdat = ft_preprocessing(cfg,data);
+
+            % clear cfg
+            clear cfg
+
+            % interpolation cfg
+            cfg.method     = 'spline';
+            cfg.badchannel = tmpdat.label(chns2interp);
+            cfg.trials     = 'all';
+            cfg.lambda     = 1e-5; 
+            cfg.order      = 4; 
+            cfg.elec       = config.elec;
+
+            % interpolate
+            tmpdat = ft_channelrepair(cfg,tmpdat);
+
+            % clear cfg
+            clear cfg chns2interp
+
+        %%  replace interpolated data
+
+            for j = 1:length(data.trial)
+
+                data.trial{j}(1:60,:) = tmpdat.trial{j};
+
+            end; clear j
+
+            % clear temporary data
+            clear tmpdat
+
+        %%  keep channel x trial artifacts  
+
+            data.ArtDect.channels_x_trials = config.ArtDect.channels_x_trials;            
+
+        %%  save data
+            data_EEG = data; clear data;
+            save([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_EEG_Rlm_Fhl_rdSeg_Art'],'data_EEG')
+
+%         check data
+%         data_EEG = THG_label_artefacts_20110916(data_EEG);
+%         figure; imagesc(config.ArtDect.channels_x_trials)
+
+            % plot resulting trials
+
+            if setting.plot == 1
+                h = figure;
+                data = NaN(64,size(data_EEG.trial,1),5540);
+                for indTrial = 1:size(data_EEG.trial,2)
+                    data(:,indTrial,1:size(data_EEG.trial{indTrial},2))= data_EEG.trial{indTrial};
+                    subplot(8,8,indTrial);
+                    imagesc(squeeze(data_EEG.trial{indTrial}));
+                end; suptitle(['processing ID ' IDs{id}, ' Run: ', num2str(iRun)]);
+                close(h);
+            end
+
+        %% create overview of percentage of excluded trials, update config.trl
+
+            tmp_trialNumbers(1,1) = length(data_EEG.trial);
+            tmp_trialNumbers(2,1) = 64; % maximum amount of trials
+            tmp_trialNumbers(3,1) = (100-((tmp_trialNumbers(1,1)/tmp_trialNumbers(2,1))*100)); % percentage of excluded trials
+            overview_trial_numbers{id, iRun} = tmp_trialNumbers;
+            config.overview_trial_numbers = tmp_trialNumbers;
+
+        %%  update config.trl with trigger codes
+
+            % add info whether trial is included in final set; 1 = trial containing no artifact
+            config.trl(:,4)                 = 0;
+            config.trl(nartfree_trials, 4)  = 1;
+            % add trial number
+            config.trl(:,5) = 0;
+            config.trl(nartfree_trials, 5)  = nartfree_trials;
+            % load marker file
+            mrk = config.mrk;
+            for i = 1:size(mrk,2)
+               mrk_val{i,1} = mrk(1,i).value;
+            end
+            % add onset and offset triggers
+            config.trl(:,6) = 17;
+            config.trl(:,7) = 64;
+            % save config
+            save([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config'],'config')
+            
+        catch ME
+            warning('Error occured. Please check.');
+            rethrow(ME)
+            fprintf('\n ID not availble! Skip! \n')
+            ID_unavailable{id,iRun} = (IDs{id});
+        end % try...catch
+    end % run
+end % id

code/08_STSW_assignConditionsToData.m

@@ -0,0 +1,129 @@
+% Add behavioral information to the EEG data.
+
+% 180123 | written by JQK
+% 180327 | including OAs, check for existing data, changed behavioral matrix
+
+clear all; close all; pack; clc;
+
+%% pathdef
+
+if ismac
+    pn.study        = '/Volumes/LNDG/Projects/StateSwitch/';
+    pn.eeg_root     = [pn.study, 'dynamic/data/eeg/task/A_preproc/SA_preproc_study/'];
+    pn.dynamic_In   = [pn.eeg_root, 'B_data/B_EEG_ET_ByRun/'];
+    pn.triggerTiming= [pn.eeg_root, 'C_figures/D_TriggerTiming/'];
+    pn.EEG          = [pn.eeg_root, 'B_data/C_EEG_FT/']; %mkdir(pn.EEG);
+    pn.History      = [pn.eeg_root, 'B_data/D_History/']; %mkdir(pn.History);
+    pn.logs         = [pn.eeg_root, 'Y_logs/H_ArtCorr/'];
+    % add ConMemEEG tools
+    pn.MWBtools     = [pn.eeg_root, 'T_tools/fnct_MWB/'];           addpath(genpath(pn.MWBtools));
+    pn.THGtools     = [pn.eeg_root, 'T_tools/fnct_THG/'];           addpath(genpath(pn.THGtools));
+    pn.commontools  = [pn.eeg_root, 'T_tools/fnct_common/'];        addpath(genpath(pn.commontools));
+    pn.fnct_JQK     = [pn.eeg_root, 'T_tools/fnct_JQK/'];           addpath(genpath(pn.fnct_JQK));
+    pn.FT           = [pn.eeg_root, 'T_tools/fieldtrip-20170904/']; addpath(pn.FT); ft_defaults;
+    pn.helper       = [pn.eeg_root, 'A_scripts/helper/'];           addpath(pn.helper);
+else
+    pn.root         = '/home/mpib/LNDG/StateSwitch/WIP_eeg/SA_preproc_study/';
+    pn.EEG          = [pn.root, 'B_data/C_EEG_FT/'];
+    pn.History      = [pn.root, 'B_data/D_History/'];
+    pn.triggerTiming= [pn.root, 'C_figures/D_TriggerTiming/'];
+    pn.logs         = [pn.root, 'Y_logs/H_ArtCorr/'];
+end
+
+%% define Condition & IDs for preprocessing
+
+condEEG = 'dynamic';
+
+%% define IDs for segmentation
+
+% N = 47 YAs + 53 OAs;
+IDs = {'1117';'1118';'1120';'1124';'1126';'1131';'1132';'1135';'1136';'1138';...
+    '1144';'1151';'1158';'1160';'1163';'1164';'1167';'1169';'1172';'1173';...
+    '1178';'1182';'1215';'1216';'1219';'1221';'1223';'1227';'1228';'1233';...
+    '1234';'1237';'1239';'1240';'1243';'1245';'1247';'1250';'1252';'1257';...
+    '1261';'1265';'1266';'1268';'1270';'1276';'1281';...
+    '2104';'2107';'2108';'2112';'2118';'2120';'2121';'2123';'2125';'2129';...
+    '2130';'2131';'2132';'2133';'2134';'2135';'2139';'2140';'2145';'2147';...
+    '2149';'2157';'2160';'2201';'2202';'2203';'2205';'2206';'2209';'2210';...
+    '2211';'2213';'2214';'2215';'2216';'2217';'2219';'2222';'2224';'2226';...
+    '2227';'2236';'2237';'2238';'2241';'2244';'2246';'2248';'2250';'2251';...
+    '2252';'2258';'2261'};
+
+%% load behavioral data
+
+behavData = load('/Volumes/LNDG/Projects/StateSwitch/dynamic/data/behavior/STSW_dynamic/A_MergeIndividualData/B_data/SS_MergedDynamic_EEG_MRI_YA_09-Mar-2018.mat');
+
+%% extract data by subject
+
+for id = 1:numel(IDs)
+    if 1%~exist([pn.EEG, IDs{id}, '_', condEEG, '_EEG_Rlm_Fhl_rdSeg_Art_Cond.mat'])
+    for iRun = 1:4
+
+        %% load preprocessed EEG data
+
+        tmp = [];
+        tmp.clock = tic; % set tic for processing duration
+
+        load([pn.History, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_config'],'config');
+        load([pn.EEG, IDs{id}, '_r',num2str(iRun), '_', condEEG, '_EEG_Rlm_Fhl_rdSeg_Art'],'data_EEG');
+
+        data = data_EEG; clear data_EEG;
+
+        %% remove extraneous channels
+
+        rem             = [];
+        rem.channel     = {'all','-IOR','-LHEOG','-RHEOG','-A1'};
+        rem.demean      = 'no';
+        
+        dataNew{iRun} = ft_preprocessing(rem,data); clear rem data;
+        
+        %% get behavioral data (incl. condition assignment)
+
+        behavIdx = find(strcmp(behavData.IDs_all, IDs{id})); % get index of current subject in behavioral matrix
+        SubjBehavData.Atts          = behavData.MergedDataEEG.Atts(:,behavIdx);
+        SubjBehavData.StateOrders   = behavData.MergedDataEEG.StateOrders(:,behavIdx);
+        SubjBehavData.RTs           = behavData.MergedDataEEG.RTs(:,behavIdx);
+        SubjBehavData.Accs          = behavData.MergedDataEEG.Accs(:,behavIdx);
+        SubjBehavData.CuedAttributes = reshape(cat(1,behavData.MergedDataEEG.expInfo{behavIdx}.AttCuesRun{:})',256,1);
+        
+        % create matrix with crucial condition info
+        CondData = NaN(256,7);
+        for indTrial = 1:256
+           tmp_cuedTrial = SubjBehavData.CuedAttributes(indTrial);
+           if ismember(1,tmp_cuedTrial{:})
+               CondData(indTrial,1) = 1;
+           end
+           if ismember(2,tmp_cuedTrial{:})
+               CondData(indTrial,2) = 1;
+           end
+           if ismember(3,tmp_cuedTrial{:})
+               CondData(indTrial,3) = 1;
+           end
+           if ismember(4,tmp_cuedTrial{:})
+               CondData(indTrial,4) = 1;
+           end
+        end
+        CondData(:,5) = repmat(1:8,1,32)'; % serial position
+        CondData(:,6) = SubjBehavData.Atts;
+        CondData(:,7) = 1:256;
+        CondData(:,8) = SubjBehavData.StateOrders;
+        CondData(:,9) = SubjBehavData.RTs;
+        CondData(:,10) = SubjBehavData.Accs;
+        % select trials from current run
+        CurTrials{iRun} = CondData(((iRun-1)*64)+1:iRun*64,:);
+        % remove trials that were not retained during preprocessing
+        removedTrials = []; removedTrials = setdiff(1:64,config.nartfree_trials);
+        CurTrials{iRun}(removedTrials,:) = [];
+        
+    end; clear CondData;
+    %% combine data across runs
+    data = ft_appenddata([], dataNew{1}, dataNew{2}, dataNew{3}, dataNew{4});
+    data.TrlInfo = cat(1,CurTrials{:});
+    data.TrlInfoLabels = {'Att1 cued', 'Att2 cued', 'Att3 cued', 'Att4 cued', ...
+        'Serial position', 'WinningAtt', 'Continuous position', 'Cue Dimensionality',...
+        'RT', 'Acc'};
+    save([pn.EEG, IDs{id}, '_', condEEG, '_EEG_Rlm_Fhl_rdSeg_Art_Cond'],'data');
+    else
+        disp(['Data for ID ',IDs{id}, ' already exists!']);
+    end; % if exist
+end