%% SCRIPT DESCRIPTION % This script combineds the results of the amicas A,B and C. It takes this information and using ICAtest it determines which components were % and were not replicated across and flags them, along with time sections with high standard deviation. Beta,theta and Alpha power is calculated % arcoss time and is marked but not primed for removal. % %% From Config key_strswap Description %----------------------------------------------------------------------------------------------------------------------------- % in_path = [in_path] Relative path to input data files assuming cd = work_path % recur_sec = [recur_sec] Recurrence (sec) for the current artifact detection epoching (e.g. 1) % limit_sec = [limit_sec] Limits (sec) for the current artifact detection epoching (e.g. [0 1]) % sd_t_meth = [sd_t_meth] Method used for flagging epochs (e.g. 'q' (quantiles), or 'na' (default)) % sd_t_vals = [sd_t_vals] Percentage trim for confidence intervals during epoch standard deviation criteria (e.g. [.3 .7]) % sd_t_o = [sd_t_o] z threshold for flagging epochs during standard deviation criteria (e.g. 6) % sd_t_f_meth = [sd_t_f_meth] Fixed method used for flagging epochs (should be 'fixed') % sd_t_f_vals = [sd_t_f_vals] Percentage trim for confidence intervals during epoch standard deviation criteria (e.g. [.3 .7], leave empty for 'fixed') % sd_t_f_o = [sd_t_f_o] z threshold for flagging epochs during fixed standard deviation criteria (e.g. .2) % min_gap_ms = [min_gap_ms] Minimum time (ms) to allow between periods marked for rejection % icr_crit = [icr_crit] False discovery rate of component replicatablilty, measure of sensitivity % % dip_vol = [dip_vol] file path/name of the standard volume file used in dipole fitting (e.g. derivatives/lossless/code/misc/standard_vol.mat). % dip_mri = [dip_mri] file path/name of the standard mri file used in dipole fitting (e.g. derivatives/lossless/code/misc/standard_mri.mat). % dip_elc = [dip_elc] file path/name of the standard electrode file used in dipole fitting (derivatives/lossless/code/misc/standard_1020.elc). % dip_transmat = [dip_transmat] Transformation matrix to apply during dipole fit (e.g. [0 0 0 0 0 0 1 1 1]) % cor_transmat = [cor_transmat] The transformation matrix to apply between dipole fitting to the standard head and performing ICLabel to the EEGLAB default orientation (e.g.[0,0,0,0,0,1.57,1,1,1]) % cor_nosedir = [cor_nosedir] Nose direction of montage to apply between dipole fitting to the standard head and performing ICLabel to the EEGLAB default orientation (e.g. +Y) % out_path = [out_path] Relative path to output data files assuming cd = work_path %% LOAD COMPART DATASET %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Loads the data file and the first amica model. logging_log('INFO', 'Loading set file: [batch_dfn,_,-1]_compart_data.set...'); tic; EEG = pop_loadset('filepath','[in_path]','filename','[batch_dfn,_,-1]_compart_data.set'); EEG = eeg_checkset( EEG ); logging_log('INFO', 'TIME TO: LOAD DATASET...'); toc % Set up the ICA channels EEG.icachansind = marks_label2index(EEG.marks.chan_info,{'manual','rank'},'indices','invert','on'); % Load Amica A logging_log('INFO', 'Loading Amica ABC...'); tic; EEG.amica(2).models=loadmodout15('[in_path]/[batch_dfn,_,-1]_amicaout_A'); EEG.icaweights=EEG.amica(2).models.W; EEG.icasphere=EEG.amica(2).models.S(1:EEG.amica(2).models.num_pcs,:); %% PASS EEG DATA THROUGH THE AMICA MATRIX %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Loads the B and C amica models. Amica A is used to fill in the ICA information. % Instead of EEG=eeg_checkset(EEG); tmpdata = EEG.data(EEG.icachansind,:); tmpindices = find(~sum(isnan(tmpdata))); EEG.icaact = zeros(size(EEG.icaweights,1), size(tmpdata,2)); EEG.icaact(:) = NaN; EEG.icaact(:,tmpindices) = (EEG.icaweights*EEG.icasphere)*tmpdata(:,tmpindices); clear tmpdata; % Create winv for topo plots manually EEG.icawinv = []; EEG.icawinv = pinv(EEG.icaweights * EEG.icasphere); % Load Amica B EEG.amica(3).models=loadmodout15('[in_path]/[batch_dfn,_,-1]_amicaout_B'); % Load Amica C EEG.amica(4).models=loadmodout15('[in_path]/[batch_dfn,_,-1]_amicaout_C'); logging_log('INFO', 'TIME TO: LOAD ICA ABC INFO...'); toc %% UPDATE REJECTION STRUCTURE... %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Creating flag marks to based on the log liklyhodd models found for each of the ICA models % Flags for amicas A B and C logging_log('INFO', 'Updating latflaginfo structure...'); tic; % A logl_A_flags=zeros(size(EEG.data(1,:))); logl_A_flag_inds=marks_label2index(EEG.marks.time_info,{'manual'},'indices','invert','on'); % Check to see if the overwrite will be successful. If the overwirte is a different size then some of the old data will be retained. if length(EEG.amica(2).models.Lt) ~= length(logl_A_flag_inds) EEG.amica(2).models.Lt = EEG.amica(2).models.Lt(1:length(logl_A_flag_inds)); warning('YOU FORGOT TO DELETE YOUR OLD AMICA FILES AND THE NEW CHANNEL SIZES DO NOT MATCH UP'); end logl_A_flags(logl_A_flag_inds)=(EEG.amica(2).models.Lt-min(EEG.amica(2).models.Lt))/(max(EEG.amica(2).models.Lt)-min(EEG.amica(2).models.Lt)); EEG.marks=marks_add_label(EEG.marks,'time_info', ... {'logl_A',[1,0.45,0],logl_A_flags}); clear logl_A_flags; % B logl_B_flags=zeros(size(EEG.data(1,:))); logl_B_flag_inds=marks_label2index(EEG.marks.time_info,{'manual'},'indices','invert','on'); logl_B_flags(logl_B_flag_inds)=(EEG.amica(3).models.Lt-min(EEG.amica(3).models.Lt))/(max(EEG.amica(3).models.Lt)-min(EEG.amica(3).models.Lt)); EEG.marks=marks_add_label(EEG.marks,'time_info', ... {'logl_B',[1,0.45,0],logl_B_flags}); clear logl_B_flags; % C logl_C_flags=zeros(size(EEG.data(1,:))); logl_C_flag_inds=marks_label2index(EEG.marks.time_info,{'manual'},'indices','invert','on'); logl_C_flags(logl_C_flag_inds)=(EEG.amica(4).models.Lt-min(EEG.amica(4).models.Lt))/(max(EEG.amica(4).models.Lt)-min(EEG.amica(4).models.Lt)); EEG.marks=marks_add_label(EEG.marks,'time_info', ... {'logl_C',[1,0.45,0],logl_C_flags}); clear logl_C_flags; logging_log('INFO', 'TIME TO: UPDATE REJECTION STRUCTURE...'); toc %% IDENTIFY RELIABLE COMPONENTS USING ISCTEST %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Gives the A variable of each of the amicas to the function icatest. Returned are the indecies of the % Components that were not replicated across all 3 of the models. tic; % Clear a in case of local models are saved already a = []; for i=2:length(EEG.amica); a(:,:,i-1)=EEG.amica(i).models(1).A; end % Calculating the replicated component flag % Replicate criteria changed from 0.05,0.05 to 0.1,0.1 to be less sensitive ca=isctest(a,[icr_crit],[icr_crit],'mixing'); if isempty(ca) mess = ['ISCTest failed to return reasonable values. ISCTest is ' ... 'meant to find duplicate channels in amica_{A,B,C} if this ' ... 'failed consider investigating the datafile']; logging_log('ERROR', mess); error(mess); end icr_inds=setdiff(1:min(size(EEG.amica(2).models(1).A)),ca(:,1)); % Creating the replicated component flag logging_log('INFO', 'Updating chflaginfo structure...'); icr_comp_flags=zeros(min(size(EEG.amica(2).models(1).A)),1); icr_comp_flags(icr_inds)=1; EEG.marks=marks_add_label(EEG.marks,'comp_info', ... {'ic_rt',[0,0,1],[0,0,1],-1,icr_comp_flags}); logging_log('INFO', 'TIME TO: EDIT FIND REPLICATED COMPONENTS AND UPDATE COMPINFO STRUCT...'); toc % Window the continuous data logging_log('INFO', 'Windowing the continous data...'); tic; EEG=marks_continuous2epochs(EEG,'recurrence',[recur_sec],'limits',[limit_sec]); logging_log('INFO', 'TIME TO: WINDOW THE CONTINUOUS DATA...'); toc %% CALCULATE ICAACT SD %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Similarly to the IC_SD flag created in compart. Marks sections of time where the components % have an abnormally high standard deviation. % Calculating the replicated component flag logging_log('INFO', 'Calculating the IC sd array for window criteria...'); tic; epoch_inds=marks_label2index(EEG.marks.time_info,{'manual'},'indexes','invert','on'); [EEG,icaact_sd2_t]=chan_variance(EEG,'data_field','icaact', ... 'epoch_inds',epoch_inds, ... 'varmeasure','absmean', ... 'detrend','off', ... 'plot_figs','off'); % Creating the ic_sd component flag logging_log('INFO', 'Assessing window icaact sd distributions...') [~,flag_sd_t_inds]=marks_array2flags(icaact_sd2_t, ... 'flag_dim','col', ... 'init_method','[sd_t_meth]', ... 'init_vals',[sd_t_vals], ... 'init_crit',[sd_t_o], ... 'flag_method','[sd_t_f_meth]', ... 'flag_vals',[ [sd_t_f_vals] ], ... 'flag_crit',[sd_t_f_o], ... 'plot_figs','off'); % Updating rejection structure logging_log('INFO', 'Updating latflaginfo structure...'); icsd_epoch_flags=zeros(size(EEG.data(1,:,:))); icsd_epoch_flags(1,:,epoch_inds(flag_sd_t_inds))=1; EEG.marks=marks_add_label(EEG.marks,'time_info', ... {'ic_sd2',[.4,.6,.6],icsd_epoch_flags}); logging_log('INFO', 'TIME TO: EDIT FIND IC_SD AND UPDATE TIMEINFO REJECTION STRUCT...'); toc clear icsd_epoch_flags; %% CALCULATE DELTA/THETA ICAACT POWER %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Calculates and flags times of high Theta/Alpha Power % Calculating the Delta/Theta component flag logging_log('INFO', 'Calculating the delta/theta IC power array for window criteria...'); tic; [EEG,icaact_dt_t]=chan_variance(EEG,'data_field','icaact', ... 'epoch_inds',epoch_inds, ... 'detrend','on', ... 'varmeasure','spect', ... 'spectrange',[0 7], ... 'plot_figs','off'); % Creating the Delta/Theta component flag logging_log('INFO', 'Assessing window icaact sd distributions...') [~,flag_dt_t_inds]=marks_array2flags(icaact_dt_t, ... 'flag_dim','col', ... 'init_method','[sd_t_meth]', ... 'init_vals',[sd_t_vals], ... 'init_crit',[sd_t_o], ... 'flag_method','[sd_t_f_meth]', ... 'flag_vals',[ [sd_t_f_vals] ], ... 'flag_crit',[sd_t_f_o], ... 'plot_figs','off'); % Updating rejection structure logging_log('INFO', 'Updating latflaginfo structure...'); icdt_epoch_flags=zeros(size(EEG.data(1,:,:))); icdt_epoch_flags(1,:,epoch_inds(flag_dt_t_inds))=1; EEG.marks=marks_add_label(EEG.marks,'time_info', ... {'ic_dt',[1,1,0],icdt_epoch_flags}); logging_log('INFO', 'TIME TO: CALCULATE DELTA/THETA ICAACT POWER AND UPDATE MARKS...'); toc clear icdt_epoch_flags; %% CALCULATE ALPHA ICAACT POWER %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Calculates and flags times of high Alpha Power % Calculating the Alpha component flag logging_log('INFO', 'Calculating the alpha IC power array for window criteria...'); tic; [EEG,icaact_a_t]=chan_variance(EEG,'data_field','icaact', ... 'epoch_inds',epoch_inds, ... 'detrend','on', ... 'varmeasure','spect', ... 'spectrange',[8 13], ... 'plot_figs','off'); % Creating the Alpha component flag logging_log('INFO', 'Assessing window icaact sd distributions...') [~,flag_a_t_inds]=marks_array2flags(icaact_a_t, ... 'flag_dim','col', ... 'init_method','[sd_t_meth]', ... 'init_vals',[sd_t_vals], ... 'init_crit',[sd_t_o], ... 'flag_method','[sd_t_f_meth]', ... 'flag_vals',[ [sd_t_f_vals] ], ... 'flag_crit',[sd_t_f_o], ... 'plot_figs','off'); % Updating rejection structure logging_log('INFO', 'Updating latflaginfo structure...'); ica_epoch_flags=zeros(size(EEG.data(1,:,:))); ica_epoch_flags(1,:,epoch_inds(flag_a_t_inds))=1; EEG.marks=marks_add_label(EEG.marks,'time_info', ... {'ic_a',[1,1,0],ica_epoch_flags}); logging_log('INFO', 'TIME TO: CALCULATE ALPHA ICAACT POWER AND UPDATE MARKS...'); toc clear ica_epoch_flags; %% CALCULATE BETA ICAACT POWER %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Calculates and flags times of high Beta Power. % Note that at the moment this mark is used for removals. % Calculating the Beta component flag logging_log('INFO', 'Calculating the beta IC power array for window criteria...'); tic; [EEG,icaact_b_t]=chan_variance(EEG,'data_field','icaact', ... 'epoch_inds',epoch_inds, ... 'detrend','on', ... 'varmeasure','spect', ... 'spectrange',[14 30], ... 'plot_figs','off'); logging_log('INFO', 'TIME TO: CALCULATE BETA ICAACT POWER...'); % Creating the Beta component flag logging_log('INFO', 'Assessing window icaact sd distributions...') [~,flag_b_t_inds]=marks_array2flags(icaact_b_t, ... 'flag_dim','col', ... 'init_method','[sd_t_meth]', ... 'init_vals',[sd_t_vals], ... 'init_crit',[sd_t_o], ... 'flag_method','[sd_t_f_meth]', ... 'flag_vals',[ [sd_t_f_vals] ], ... 'flag_crit',[sd_t_f_o], ... 'plot_figs','off'); logging_log('INFO', 'TIME TO: CREATE WINDOW CRITERIA VECTOR...'); % Updating rejection structure logging_log('INFO', 'Updating latflaginfo structure...'); icb_epoch_flags=zeros(size(EEG.data(1,:,:))); icb_epoch_flags(1,:,epoch_inds(flag_b_t_inds))=1; EEG.marks=marks_add_label(EEG.marks,'time_info', ... {'ic_b',[0,1,1],icb_epoch_flags}); logging_log('INFO', 'TIME TO: CALCULATE BETA ICAACT POWER AND UPDATE MARKS...'); toc clear icb_epoch_flags; %% CALCULATE LOW GAMMA ICAACT POWER %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Calculates and flags times of high Low-Gamma Power. % Note that at the moment this mark is used for removals. % Calculating the Beta component flag logging_log('INFO', 'Calculating the low-gamma IC power array for window criteria...'); tic; [EEG,icaact_lg_t]=chan_variance(EEG,'data_field','icaact', ... 'epoch_inds',epoch_inds, ... 'detrend','on', ... 'varmeasure','spect', ... 'spectrange',[31 60], ... 'plot_figs','off'); logging_log('INFO', 'TIME TO: CALCULATE LOW-GAMMA ICAACT POWER...'); % Creating the LOW-GAMMA component flag logging_log('INFO', 'Assessing window icaact sd distributions...') [~,flag_lg_t_inds]=marks_array2flags(icaact_lg_t, ... 'flag_dim','col', ... 'init_method','[sd_t_meth]', ... 'init_vals',[sd_t_vals], ... 'init_crit',[sd_t_o], ... 'flag_method','[sd_t_f_meth]', ... 'flag_vals',[ [sd_t_f_vals] ], ... 'flag_crit',[sd_t_f_o], ... 'plot_figs','off'); logging_log('INFO', 'TIME TO: CREATE WINDOW CRITERIA VECTOR...'); % Updating rejection structure logging_log('INFO', 'Updating latflaginfo structure...'); iclg_epoch_flags=zeros(size(EEG.data(1,:,:))); iclg_epoch_flags(1,:,epoch_inds(flag_lg_t_inds))=1; EEG.marks=marks_add_label(EEG.marks,'time_info', ... {'ic_lg',[0,1,1],iclg_epoch_flags}); logging_log('INFO', 'TIME TO: CALCULATE LOW-GAMMA ICAACT POWER AND UPDATE MARKS...'); toc clear iclg_epoch_flags; %% CALCULATE HIGH GAMMA ICAACT POWER %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Calculates and flags times of high High-Gamma Power. % Note that at the moment this mark is used for removals. % Calculating the Beta component flag logging_log('INFO', 'Calculating the high-gamma IC power array for window criteria...'); tic; [EEG,icaact_hg_t]=chan_variance(EEG,'data_field','icaact', ... 'epoch_inds',epoch_inds, ... 'detrend','on', ... 'varmeasure','spect', ... 'spectrange',[61 100], ... 'plot_figs','off'); logging_log('INFO', 'TIME TO: CALCULATE HIGH-GAMMA ICAACT POWER...'); % Creating the HIGH-GAMMA component flag logging_log('INFO', 'Assessing window icaact sd distributions...') [~,flag_hg_t_inds]=marks_array2flags(icaact_hg_t, ... 'flag_dim','col', ... 'init_method','[sd_t_meth]', ... 'init_vals',[sd_t_vals], ... 'init_crit',[sd_t_o], ... 'flag_method','[sd_t_f_meth]', ... 'flag_vals',[ [sd_t_f_vals] ], ... 'flag_crit',[sd_t_f_o], ... 'plot_figs','off'); logging_log('INFO', 'TIME TO: CREATE WINDOW CRITERIA VECTOR...'); % Updating rejection structure logging_log('INFO', 'Updating latflaginfo structure...'); ichg_epoch_flags=zeros(size(EEG.data(1,:,:))); ichg_epoch_flags(1,:,epoch_inds(flag_hg_t_inds))=1; EEG.marks=marks_add_label(EEG.marks,'time_info', ... {'ic_hg',[0,1,1],ichg_epoch_flags}); logging_log('INFO', 'TIME TO: CALCULATE HIGH-GAMMA ICAACT POWER AND UPDATE MARKS...'); toc clear ichg_epoch_flags; %% MARKS UPDATE %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Merges the epochs back to continuous data. Flags marks gaps and combines time % the new time marks into the manual marks. % Concat epochs back into continuous data logging_log('INFO', 'Concatenating windowed data...'); tic; EEG=marks_epochs2continuous(EEG); EEG=eeg_checkset(EEG,'eventconsistency'); logging_log('INFO', 'TIME TO: CONCATENATE THE WINDOWED DATA INTO CONTINUOUS DATA...'); toc % Mark flag gaps between time marks tic; EEG=pop_marks_flag_gap(EEG,{'manual','mark_gap','ic_sd2'}, ... [min_gap_ms],'mark_gap',[.8,.8,.8],'offsets',[0 0],'ref_point','both'); logging_log('INFO', 'TIME TO: MARK FLAG GAPS THAT ARE < 2 SECONDS...'); toc % Combine the removal marks into manual tic; EEG=pop_marks_merge_labels(EEG,'time_info',{'manual','ic_sd2','mark_gap'},'target_label','manual'); logging_log('INFO', 'TIME TO: COMBINE MARKS STRUCTURE INTO MANUAL FLAGS...'); toc % begin prep for dipole fit logging_log('INFO','Checking for JVM...'); javachk('jvm'); % TODO(brad) this isn't added for some reason addpath([pwd '/derivatives/lossless/code/dependencies/eeglab_asr_amica/plugins/Fieldtrip-lite141209/src']) %% DIPFIT %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Does the dipole fitting procedures. This takes a substantial amount of % time. % TODO: Make sure this is all needed, seems to be repeating things logging_log('INFO','FITTING DIPOLES...'); tic; EEG = pop_dipfit_settings( EEG, 'hdmfile','[dip_vol]', ... 'coordformat','MNI', ... 'mrifile','[dip_mri]', ... 'chanfile','[dip_elc]', ... 'coord_transform',[dip_transmat], ... 'chansel',[1:EEG.nbchan] ); % OLD SINGLE THREADED SETUP % EEG = pop_multifit(EEG, [1:EEG.nbchan] ,'threshold',100,'plotopt',{'normlen' 'on'}); % NEW MULTI-THREADED SETUP EEG = pop_multifit_bucanl_parallel(EEG, [1:EEG.nbchan], 8,'threshold',100,'plotopt',{'normlen' 'on'}); EEG=eeg_checkset(EEG); logging_log('INFO','TIME TO: FIT DIPOLES...'); toc %% MARK ICs FOR REJECTION BASED ON RV %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Creates marks based on dipfit results. Dipfit marks are added % for rejection while ICMARC are just for reference. % IC RESIDUAL VARIANCE FLAG %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Residuals greater than 15% are flagged - can be done here as it only references % dipoles and their values. %logging_log('INFO','UPDATE MARKS STRUCTURE WITH IC_RV...'); %tic; % Dipole criteria changed from 0.15 to 0.2 to be less sensitive %EEG.reject.gcompreject(find([EEG.dipfit.model.rv]>.15))=1; %EEG.reject.gcompreject(find(isnan([EEG.dipfit.model.rv])))=1; %EEG.marks=marks_add_label(EEG.marks,'comp_info', ... % {'ic_rv',[1,0,0],[1,0,0],-1,EEG.reject.gcompreject'}); %logging_log('INFO','TIME TO: UPDATE MARKS STRUCTURE...'); %toc % COMBINE MARKS STRUCTURE INTO MANUAL FLAGS... %tic; %EEG=pop_marks_merge_labels(EEG,'comp_info',{'manual','ic_rv'},'target_label','manual'); %logging_log('INFO','TIME TO: COMBINE MARKS STRUCTURE INTO MANUAL FLAGS...'); %toc % Adjust montage %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % This will pretty much always be dataset specific - so be careful EEG = warp_locs( EEG, '[dip_elc]', ... 'transform',[cor_transmat], ... 'manual','off'); EEG.chaninfo.nosedir='[cor_nosedir]'; EEG=eeg_checkset(EEG); %START TEMPORARY MOVE TO QC SCRIP % %% Use the ICLabel plugin % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % ICLabel is a database style plugin that classifies ICA components into the % % following seven categories % % % (Brain)(Muscle)(Eye)(Heart)(Line Noise)(Channel Noise)(Other) % % ( 1 )( 2 )( 3 )( 4 )( 5 )( 6 )( 7 ) % % % Colour structure and groupings are as follows: % % ---------- % % Brain -> [0, 0.4, 0] % % ---------- % % Muscle -> [0.4, 0, 0] % % Eye -> [0.8, 0, 0] % % Heart -> [0.8, 0.4, 0] % % ---------- % % Line Noise -> [0.2, 0.4, 0.8] % % Channel noise -> [0.4, 0, 0.8] % % ---------- % % Other -> [0.8, 0.8, 0] % % ---------- % % % Create lookup table of the above information. % pairs = {'brain',[0, 0.4, 0]; % 'muscle',[0.4,0,0]; % 'eye',[0.8,0,0]; % 'heart',[0.8,0.4,0]; % 'line_noise',[0.2,0.4,0.8]; % 'chan_noise',[0.4,0,0.8]; % 'other',[0.8,0.8,0]} % lookupTable = cell2struct(pairs, {'name', 'colour'}, 2); % % % Below code is how to call proper label lookup % % [prob, classind] = max(EEG.etc.ic_classification.ICLabel.classifications(1, :)); % % Below line should give you number of components % % length(EEG.etc.ic_classification.ICLabel.classifications); % % % Run ICLabel % EEG = iclabel(EEG); % % % Extract labels from the structure. % rawLabels = zeros([1,length(EEG.etc.ic_classification.ICLabel.classifications)]); % for index = 1:length(EEG.etc.ic_classification.ICLabel.classifications) % [p,i] = max(EEG.etc.ic_classification.ICLabel.classifications(index, :)); % rawLabels(index) = i; % end % % %ICLabel Exporting to marks % %%%%%%%%%%%%%%%%%%%%%%%%%%%% % logging_log('INFO','UPDATE MARKS STRUCTURE WITH ICLABEL BRAIN COMPS...'); % tic; % for compStage = 1:7 % ICLabel_pass= zeros([1,length(rawLabels)]); % ICLabel_pass(find(rawLabels==compStage)) = 1; % ICLabel_pass = ICLabel_pass'; % EEG.marks=marks_add_label(EEG.marks,'comp_info', ... % {lookupTable(compStage).name,lookupTable(compStage).colour, ... % lookupTable(compStage).colour,1,ICLabel_pass}); % end % logging_log('INFO','TIME TO: UPDATE MARKS STRUCTURE...'); % toc %END TEMPORARY MOVE TO QC SCRIP %% SAVE LOSSLESS DERIVATIVE *_eeg_ll.set FILE (same name as the input file from the % BIDS root directory [but located in the derivatives/lossless directory]). %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% logging_log('INFO','Saving set file: [out_path]/[batch_dfn,_,-1]_eeg_ll.set...'); tic; EEG = pop_saveset( EEG, 'filename','[out_path]/[batch_dfn,_,-1]_eeg_ll.set'); logging_log('INFO','TIME TO: SAVE DIP DATA SET...'); toc % Save diagnostic arrays try OCTAVE_VERSION; save('-mat7-binary','[out_path]/[batch_dfn,_,-1]_icaact_sd2_t.mat','icaact_sd2_t'); catch save('[out_path]/[batch_dfn,_,-1]_icaact_sd2_t.mat','icaact_sd2_t'); end logging_log('INFO', 'Scheduler: [scheduler]'); print_chan_sample(EEG);