sets2TFR

 add other fields to tfrGND like expdesc

0001 function tfrGND=sets2TFR(set_fnames,time_window_length,time_step,use_bins,foi)
0002 
0003 % add other fields to tfrGND like expdesc
0004 n_subs=length(set_fnames);
0005 T=time_window_length; %length of time window ?? make an option
0006 
0007 % set_fnames=cell(1,n_subs);
0008 % for sub=1:n_subs,
0009 %     set_fnames{sub}=['/Users/dgroppe/Documents/MATLAB/DATA_EYES/THREESEC/eyes0' int2str(sub) '_short.set'];
0010 % end
0011 
0012 for s=1:n_subs,
0013     if s==1,
0014         global EEG; %necessary for removing artifacts
0015         %first subject, initialize everything
0016         EEG=pop_loadset(set_fnames{1});
0017         n_bins=length(EEG.bindesc);
0018         if isempty(use_bins),
0019             use_bins=1:n_bins;
0020             n_use_bins=n_bins;
0021         elseif max(use_bins)>n_bins,
0022             error('This file only has %d bins, but you requested importing bin %d.\n', ...
0023                 n_bins,max(use_bins));
0024         else
0025             n_use_bins=length(use_bins);
0026         end
0027         
0028         [n_chans, n_tpts, n_trials]=size(EEG.data);
0029         srate=EEG.srate;
0030         orig_times=EEG.times;
0031         orig_chanlocs=EEG.chanlocs;
0032         orig_bindesc=EEG.bindesc;
0033         
0034         %intialize tfrGND variable
0035         tfrGND.bindesc=cell(1,n_use_bins);
0036         for b=1:n_use_bins,
0037             tfrGND.bindesc{b}=EEG.bindesc{use_bins(b)};
0038         end
0039         tfrGND.set_fnames=set_fnames;
0040         
0041         subXbin.ftrip=cell(n_subs,n_use_bins);
0042         
0043         remove_artifact_ics(); %ICA artifact correction
0044         %in_bin=bin_membership(EEG); % figure out which epochs belong in which bins
0045         
0046         dataPP=eeglab2fieldtrip(EEG,'preprocessing','none');
0047         
0048         cfg_pre=[];
0049         % ?? make channel selection an option
0050         %cfg_pre.channel={'all','-lle','-lhz','-rle','-rhz','-A2'}; %note this only removes channel data, not their info in dataPP2.elec
0051         cfg_pre.channel={'all','-A2'}; %make this an option!!!
0052         cfg_pre.demean='yes'; %or 'no'
0053         cfg_pre.baselinewindow ='all'; % or [begin end]
0054         cfg_pre.trials = 'all'; %?? make selecting a subset of bins possible
0055         
0056         tpt0_id=find_tpt(EEG.times,0);
0057         
0058         %We have to add a trial configuration to make sure fieldtrip
0059         %doesn't think the data are really continuous.  We don't use it for
0060         %anything though (currently).
0061         for a=1:n_trials,
0062             dataPP.cfg.trl(a,:)=[1 n_tpts tpt0_id a]; %"a" just tells you what EEG.epoch the trial corresponds to. It isn't used for anything now.
0063             %get rid of possible time numerical error; milliseconds is the
0064             %lowest resolution we should ever need
0065             dataPP.time{a}=round(dataPP.time{a}*1000)/1000;
0066         end
0067         
0068         dataPP=ft_preprocessing(cfg_pre,dataPP);
0069     else
0070         EEG=pop_loadset(set_fnames{s});
0071         n_bins2=length(EEG.bindesc);
0072         [n_chans2, n_tpts2, n_trials2]=size(EEG.data);
0073         %get rid of possible time numerical error; milliseconds is the
0074         %lowest resolution we should ever need
0075         EEG.times=round(EEG.times);
0076         
0077         
0078  
0079         
0080         %Check for compatibility with previously loaded set files
0081         if EEG.srate~=srate,
0082             error('Sampling rate of file %s does not equal that of %s.\n',set_fnames{s},set_fnames{1});
0083         end
0084         if n_bins2~=n_bins,
0085             error('File %s has a different number of bins than %s.\n',set_fnames{s},set_fnames{1});
0086         end
0087         if ~isequal(EEG.bindesc,orig_bindesc)
0088             error('File %s has a different bin descriptors than %s.\n',set_fnames{s},set_fnames{1});
0089         end
0090         if n_chans2~=n_chans,
0091             error('File %s has a different number of channels (i.e., electrodes) than %s.\n',set_fnames{s},set_fnames{1});
0092         end
0093         if ~isequal(EEG.chanlocs,orig_chanlocs),
0094             error('File %s has different channel location information than %s.\n',set_fnames{s},set_fnames{1});
0095         end
0096         if n_tpts2~=n_tpts,
0097             error('File %s has a different number of time points than %s.\n',set_fnames{s},set_fnames{1});
0098         end
0099         if ~isequal(EEG.times,orig_times),
0100             error('File %s has a different time point values than %s.\n',set_fnames{s},set_fnames{1});
0101         end
0102         
0103         remove_artifact_ics(); %ICA artifact correction
0104         
0105         dataPP=eeglab2fieldtrip(EEG,'preprocessing','none');
0106         
0107         tpt0_id=find_tpt(EEG.times,0);
0108         %We have to add a trial configuration to make sure fieldtrip
0109         %doesn't think the data are really continuous.  We don't use it for
0110         %anything though (currently).
0111         for a=1:n_trials2,
0112             dataPP.cfg.trl(a,:)=[1 n_tpts tpt0_id a]; %"a" just tells you what EEG.epoch the trial corresponds to. It isn't used for anything now.
0113             %get rid of possible time numerical error; milliseconds is the
0114             %lowest resolution we should ever need
0115             dataPP.time{a}=round(dataPP.time{a}*1000)/1000;
0116         end
0117   
0118         dataPP=ft_preprocessing(cfg_pre,dataPP);
0119     end
0120     
0121     %% compute spectrogram for each bin,
0122     in_bin=bin_membership(EEG); % figure out which epochs belong in which bins
0123     
0124     %fieldtrip statistical analyses and grandaveraging expect different
0125     %conditions to be stored in different subXbin variables
0126     cfg              = [];
0127     cfg.output       = 'pow';
0128     cfg.channel      = 'all';
0129     cfg.method       = 'mtmconvol';
0130 
0131     %cfg.pad=0; %make an option ??
0132     %cfg.taper        = 'dpss'; ??
0133     cfg.taper        = 'hanning';
0134     %cfg.tapsmofrq    = 1; % multitaper option fix ??
0135   
0136     %cfg.foilim        = [0 EEG.srate/2]; %?? make an option
0137     freq_step=1;
0138     if isempty(foi),
0139         cfg.foi          = 0:freq_step:srate/2; % Spectral power will be computed at these frequencies
0140     else
0141         cfg.foi          = foi(1):freq_step:foi(end); % Spectral power will be computed at these frequencies
0142     end
0143     
0144     cfg.pad = ceil((n_tpts/srate)*freq_step)/freq_step; %"pad" is in units of seconds;
0145     %fieldtrip requires "pad" to be at least as long as the length of each
0146     %epoch of data.  Here we round up such that the pad length is a multiple
0147     %of our desired frequency step size
0148     
0149     cfg.t_ftimwin    = ones(length(cfg.foi),1).*T;   % length of time window for each frequency
0150 
0151     %find earliest possible time point to center the moving window
0152     start_time=dataPP.time{1}(1)+T/2; 
0153     start_time=ceil(start_time*srate)/srate; %you need to round in case start_time falls in between samples
0154     
0155     %latest possible time point to center the moving window
0156     end_time=dataPP.time{1}(end)-T/2; 
0157     end_time=floor(end_time*srate)/srate; %you need to round in case start_time falls in between samples
0158     cfg.toi          = start_time:time_step:end_time; %make an option ??
0159 
0160     %perform the time frequency analysis for each bin for this subject
0161     for b=1:n_use_bins,
0162       cfg.trials=find(in_bin(:,use_bins(b)));
0163       subXbin.ftrip{s,b}=ft_freqanalysis(cfg, dataPP);
0164     end
0165     
0166 end
0167 
0168 %% construct grand average with individual bins
0169 for bin=1:n_use_bins,
0170     cfg=[];
0171     cfg.keepindividual = 'yes';
0172     cmnd=['tfrGND.ftrip{bin}=ft_freqgrandaverage(cfg'];
0173     for s=1:n_subs,
0174         cmnd=[cmnd  ',subXbin.ftrip{' int2str(s) ',' int2str(bin) '}'];
0175     end
0176     cmnd=[cmnd ');'];
0177     warning('off','all');
0178     eval(cmnd);
0179     warning('on','all');
0180     tfrGND.ftrip{bin}.elec=dataPP.elec; %add electrode info, fieldtrip gets rid of it
0181 end
0182 
0183 
0184 %% construct grand average withOUT individual bins
0185 for bin=1:n_use_bins,
0186     tfrGND.gnd_powspctrm{bin}=squeeze(mean(tfrGND.ftrip{bin}.powspctrm));
0187     %Code below does this using fieldtrip
0188     %     cfg=[];
0189     %     cfg.keepindividual = 'no';
0190     %     cmnd=['grands{bin}=ft_freqgrandaverage(cfg'];
0191     %     for s=1:n_subs,
0192     %         cmnd=[cmnd  ',subXbin.ftrip{' int2str(s) ',' int2str(bin) '}'];
0193     %     end
0194     %     cmnd=[cmnd ');'];
0195     %     eval(cmnd);
0196 end
0197 tfrGND.gnd_dimord='chan_freq_time';
0198 
0199 % stuff i could add to tfrGND struct ??
0200 %     * exp_desc:<-the name of the experiment
0201 %     * filename:<-the filename of the GND variable
0202 %     * filepath:<-the filepath of the GND variable
0203 %     * saved: <-whether or not the GND variable has been saved since it
0204 %     was last modified
0205 %     * bsln_wind: <-the time window used to baseline the ERPs (in
0206 %     milliseconds)
0207 %     * indiv_subnames:<-The codenames of each subject
0208 %     * indiv_traits: <-A field where you could store traits of individual subjects (e.g., age, working memory span)
0209 %     * indiv_bin_ct: <-The number of trials each participant contributed to each bin.
0210 %     * indiv_bin_raw_ct: <-A field specific to Kutaslab data that can be
0211 %     ignored by non-Kutaslab users.
0212 %     * indiv_art_ics: <-The numeric indices of independent components
0213 %     removed from each set file before using them to compute ERPs
0214 %     * history: <-A record of commands applied to the GND variable
0215 %     * t_tests: <-The results of mass univariate t-test analyses applied to this GND variable
0216 
0217 
0218 function in_bin=bin_membership(EEG)
0219 %function in_bin=bin_membership(EEG)
0220 
0221 n_bin=length(EEG.bindesc);
0222 n_epoch=length(EEG.epoch);
0223 
0224 in_bin=zeros(n_epoch,n_bin);
0225 
0226 for a=1:n_bin,
0227     bin_name=['bin' int2str(a)];
0228     for ep=1:n_epoch,
0229         if ismember(bin_name,EEG.epoch(ep).eventtype)
0230             in_bin(ep,a)=1;
0231         end
0232     end
0233 end
0234 
0235