Component analysis
Another cleaning method which is good for HB and also for blinks is ICA. pca may also work for large artifacts such as HB. first let's read the visual data.
Contents
Looking at the data without cleaning (just filter)
cd somsens fileName='c,rfhp0.1Hz'; % we read an uncleaned file so we can see HB as a component. cfg=[]; cfg.dataset=fileName; cfg.trialdef.eventtype='TRIGGER'; cfg.trialdef.prestim=0.1; cfg.trialdef.poststim=0.5; cfg.trialdef.offset=-0.1; cfg.trialdef.visualtrig= 'visafter'; cfg.trialfun='BIUtrialfun'; cfg.trialdef.eventvalue= [222 230 240 250]; %left index finger cfg1=ft_definetrial(cfg); cfg1.demean='yes'; cfg1.baselinewindow=[-0.1,0]; cfg1.bpfilter='yes'; cfg1.bpfreq=[1 40]; cfg1.channel='MEG'; cfg1.feedback='no'; dataorig=ft_preprocessing(cfg1); % averaging cfg=[]; cfg.feedback='no'; allvis=ft_timelockanalysis(cfg,dataorig); plot(allvis.time,allvis.avg); cfgp=[]; cfgp.layout='4D248.lay'; cfgp.interactive='yes'; cfgp.zlim=[-10^-13 10^-13]; cfgp.xlim=[0.12 0.12]; figure; ft_topoplotER(cfgp,allvis); figure; cfgp.xlim=[0.19 0.19]; ft_topoplotER(cfgp,allvis); cfgp=rmfield(cfgp,'xlim'); fig1=figure; set(fig1,'Position',[0,0,800,800]); ft_multiplotER(cfgp,allvis);
evaluating trialfunction 'BIUtrialfun'
reading header...
applying digital weights in the gradiometer balancing matrix
undoing the Supine balancing
reading events from file...
applying digital weights in the gradiometer balancing matrix
undoing the Supine balancing
found 894 events
created 160 trials
the call to "ft_definetrial" took 3 seconds and an estimated 0 MB
applying digital weights in the gradiometer balancing matrix
undoing the Supine balancing
processing channel { 'A22' 'A2' 'A104' 'A241' 'A138' 'A214' 'A71' 'A26' ...
the call to "ft_preprocessing" took 7 seconds and an estimated 0 MB
the input is raw data with 248 channels and 160 trials
the call to "ft_timelockanalysis" took 0 seconds and an estimated 0 MB
reading layout from file 4D248.lay
the call to "ft_prepare_layout" took 0 seconds and an estimated 0 MB
Warning: Some points fall outside the outline, please consider using another layout
the call to "ft_topoplotTFR" took 0 seconds and an estimated 0 MB
the call to "ft_topoplotER" took 0 seconds and an estimated 0 MB
reading layout from file 4D248.lay
the call to "ft_prepare_layout" took 0 seconds and an estimated 0 MB
the call to "ft_topoplotTFR" took 0 seconds and an estimated 0 MB
the call to "ft_topoplotER" took 0 seconds and an estimated 0 MB
selection avg along dimension 1
selection dof along dimension 1
selection var along dimension 1
reading layout from file 4D248.lay
the call to "ft_prepare_layout" took 0 seconds and an estimated 0 MB
the call to "ft_multiplotER" took 1 seconds and an estimated 0 MB
Component analysis
ICA takes a lot of time. pca is faster, but!
cfgc = []; cfgc.method='pca'; cfgc.numcomponent=20; comp = ft_componentanalysis(cfgc, dataorig); %see the components and find the artifact cfgb=[]; cfgb.layout='4D248.lay'; cfgb.channel = {comp.label{1:5}}; cfgbo=ft_databrowser(cfgb,comp); % see the ica components. % look for heartbeat. it repeats almost every trial. remember it's number. % there are two ways to reject a component. reconstracting the data from % the good parts and reducing the bad component from the original data. % here is the first way of doing it. cfgrc = []; cfgrc.component = 1; % change cfgrc.feedback='no'; dataca = ft_rejectcomponent(cfgrc, comp); % here is the second % dataca = ft_rejectcomponent(cfgrc, comp, dataorig); % the second way is usefull when you only create some 10 or 20 components rather than 248 % with cfgc.numcomponent=20; % because you cannot reconstruct the data well with 20 components. % it is also usefull when you calculate the components with one piece of % data (where you have lots of blinks) in order to clean another piece of % data.
the input is raw data with 248 channels and 160 trials selecting 248 channels baseline correcting data scaling data with 1 over 0.000000 concatenating data................................................................................................................................................................ concatenated data matrix size 248x97760 starting decomposition using pca applying the mixing matrix to the sensor description the call to "ft_componentanalysis" took 1 seconds and an estimated 185 MB reading layout from file 4D248.lay the call to "ft_prepare_layout" took 1 seconds and an estimated 185 MB the input is component data with 20 components and 248 original channels detected 0 visual artifacts redrawing with viewmode component fetching data... done fetching artifacts... done preprocessing data... done the call to "ft_prepare_layout" took 0 seconds and an estimated 0 MB plotting artifacts... plotting events... plotting data... the call to "ft_prepare_layout" took 0 seconds and an estimated 0 MB plotting component topographies... done the call to "ft_databrowser" took 150 seconds and an estimated 185 MB removing 1 components keeping 19 components the call to "ft_rejectcomponent" took 1 seconds and an estimated 0 MB![]()
Component analysis on raw data
we have eye movement block in the end of the experiment
trig=readTrig_BIU('c,rfhp0.1Hz'); trig=clearTrig(trig); startt=find(trig==50,1)/1017.25; endt=find(trig==52,1)/1017.25; cfg=[]; cfg.dataset='c,rfhp0.1Hz'; cfg.trialdef.beginning=startt; cfg.trialdef.end=endt; cfg.trialfun='trialfun_raw'; % the other usefull trialfun we have are trialfun_beg and trialfun_BIU cfg1=ft_definetrial(cfg); cfg1.demean='yes';% old version was: cfg1.blc='yes'; cfg1.baselinewindow=[-0.1,0]; cfg1.lpfilter='yes'; cfg1.lpfreq=40; cfg1.channel='MEG'; MOG=ft_preprocessing(cfg1); % lets view the raw data for one channel cfgb=[]; cfgb.layout='4D248.lay'; cfgb.continuous='yes'; cfgb.event.type=''; cfgb.event.sample=1; cfgb.blocksize=3; cfgb.channel='A248'; comppic=ft_databrowser(cfgb,MOG);
Warning: 50Hz cleaning with cleanMEG pack will not be possible using the new trigger
evaluating trialfunction 'trialfun_raw'
reading header...
Warning: READ_HEADER is only a compatibility wrapper, which will soon be removed. Please instead
call FT_READ_HEADER.
applying digital weights in the gradiometer balancing matrix
undoing the Supine balancing
found 1 events
created 1 trials
the call to "ft_definetrial" took 1 seconds and an estimated 0 MB
applying digital weights in the gradiometer balancing matrix
undoing the Supine balancing
processing channel { 'A22' 'A2' 'A104' 'A241' 'A138' 'A214' 'A71' 'A26' 'A93' 'A39' 'A125' 'A20' 'A65' 'A9' 'A8' 'A95' 'A114' 'A175' 'A16' 'A228' 'A35' 'A191' 'A37' 'A170' 'A207' 'A112' 'A224' 'A82' 'A238' 'A202' 'A220' 'A28' 'A239' 'A13' 'A165' 'A204' 'A233' 'A98' 'A25' 'A70' 'A72' 'A11' 'A47' 'A160' 'A64' 'A3' 'A177' 'A63' 'A155' 'A10' 'A127' 'A67' 'A115' 'A247' 'A174' 'A194' 'A5' 'A242' 'A176' 'A78' 'A168' 'A31' 'A223' 'A245' 'A219' 'A12' 'A186' 'A105' 'A222' 'A76' 'A50' 'A188' 'A231' 'A45' 'A180' 'A99' 'A234' 'A215' 'A235' 'A181' 'A38' 'A230' 'A91' 'A212' 'A24' 'A66' 'A42' 'A96' 'A57' 'A86' 'A56' 'A116' 'A151' 'A141' 'A120' 'A189' 'A80' 'A210' 'A143' 'A113' 'A27' 'A137' 'A135' 'A167' 'A75' 'A240' 'A206' 'A107' 'A130' 'A100' 'A43' 'A200' 'A102' 'A132' 'A183' 'A199' 'A122' 'A19' 'A62' 'A21' 'A229' 'A84' 'A213' 'A55' 'A32' 'A85' 'A146' 'A58' 'A60' 'A88' 'A79' 'A169' 'A54' 'A203' 'A145' 'A103' 'A163' 'A139' 'A49' 'A166' 'A156' 'A128' 'A68' 'A159' 'A236' 'A161' 'A121' 'A4' 'A61' 'A6' 'A126' 'A14' 'A94' 'A15' 'A193' 'A150' 'A227' 'A59' 'A36' 'A225' 'A195' 'A30' 'A109' 'A172' 'A108' 'A81' 'A171' 'A218' 'A173' 'A201' 'A74' 'A29' 'A164' 'A205' 'A232' 'A69' 'A157' 'A97' 'A217' 'A101' 'A124' 'A40' 'A123' 'A153' 'A178' 'A1' 'A179' 'A33' 'A147' 'A117' 'A148' 'A87' 'A89' 'A243' 'A119' 'A52' 'A142' 'A211' 'A190' 'A53' 'A192' 'A73' 'A226' 'A136' 'A184' 'A51' 'A237' 'A77' 'A129' 'A131' 'A198' 'A197' 'A182' 'A46' 'A92' 'A41' 'A90' 'A7' 'A23' 'A83' 'A154' 'A34' 'A17' 'A18' 'A248' 'A149' 'A118' 'A208' 'A152' 'A140' 'A144' 'A209' 'A110' 'A111' 'A244' 'A185' 'A246' 'A162' 'A106' 'A187' 'A48' 'A221' 'A196' 'A133' 'A158' 'A44' 'A134' 'A216' }
reading and preprocessing
reading and preprocessing trial 1 from 1
the call to "ft_preprocessing" took 1 seconds and an estimated 0 MB
the input is raw data with 248 channels and 1 trials
detected 0 visual artifacts
redrawing with viewmode butterfly
fetching data... done
fetching artifacts... done
preprocessing data... done
plotting artifacts...
plotting events...
plotting data...
done
the call to "ft_databrowser" took 5 seconds and an estimated 0 MB
PCA
cfgp=[]; cfgc.method='pca'; compMOG = ft_componentanalysis(cfgc, MOG); % see the components and find the HB and MOG artifact % remember the numbers of the bad components and close the data browser cfgb=[]; cfgb.layout='4D248.lay'; cfgb.channel = {comp.label{1:5}}; cfgb.continuous='yes'; cfgb.event.type=''; cfgb.event.sample=1; cfgb.blocksize=3; comppic=ft_databrowser(cfgb,compMOG); % set the bad comps as the value for cfgrc.component. cfgrc = []; cfgrc.component = [1 2]; % change cfgrc.feedback='no'; dataca = ft_rejectcomponent(cfgrc, compMOG,dataorig);
the input is raw data with 248 channels and 1 trials selecting 248 channels baseline correcting data scaling data with 1 over 0.000000 concatenating data. concatenated data matrix size 248x9050 starting decomposition using pca applying the mixing matrix to the sensor description the call to "ft_componentanalysis" took 0 seconds and an estimated 0 MB reading layout from file 4D248.lay the call to "ft_prepare_layout" took 1 seconds and an estimated 0 MB the input is component data with 20 components and 248 original channels detected 0 visual artifacts redrawing with viewmode component fetching data... done fetching artifacts... done preprocessing data... done the call to "ft_prepare_layout" took 0 seconds and an estimated 0 MB plotting artifacts... plotting events... plotting data... the call to "ft_prepare_layout" took 0 seconds and an estimated 0 MB plotting component topographies... done the call to "ft_databrowser" took 15 seconds and an estimated 0 MB baseline correcting data removing 2 components keeping 18 components the call to "ft_rejectcomponent" took 4 seconds and an estimated 0 MB
reject visual trial by trial
cfg=[]; cfg.method='trial'; %trial cfg.channel='MEG'; cfg.alim=1e-12; datacln=ft_rejectvisual(cfg, dataca);
the input is raw data with 248 channels and 160 trials showing the data per trial, all channels at once trial 1 marked as GOOD trial 1 marked as GOOD 160 trials marked as GOOD, 0 trials marked as BAD 248 channels marked as GOOD, 0 channels marked as BAD no trials were removed the call to "ft_rejectvisual" took 4 seconds and an estimated 0 MB
reject visual by variance
cfg=[]; cfg.method='summary'; %trial cfg.channel='MEG'; cfg.alim=1e-12; datacln=ft_rejectvisual(cfg, dataca); cfg=[]; cfg.feedback='no'; clnAvg=ft_timelockanalysis(cfg,datacln); cfgmp=[]; cfgmp.layout='4D248.lay'; cfgmp.interactive='yes'; fig2=figure; set(fig2,'Position',[0,0,800,800]); ft_multiplotER(cfgmp,clnAvg,allvis); title ('Before (blue) and after cleaning(red)')
the input is raw data with 248 channels and 160 trials showing a summary of the data for all channels and trials computing metric [--------------------------------------------------------\] 160 trials marked as GOOD, 0 trials marked as BAD 248 channels marked as GOOD, 0 channels marked as BAD no trials were removed the call to "ft_rejectvisual" took 3 seconds and an estimated 0 MB the input is raw data with 248 channels and 160 trials the call to "ft_timelockanalysis" took 0 seconds and an estimated 0 MB selection avg along dimension 1 selection dof along dimension 1 selection var along dimension 1 selection avg along dimension 1 selection dof along dimension 1 selection var along dimension 1 reading layout from file 4D248.lay the call to "ft_prepare_layout" took 0 seconds and an estimated 0 MB the call to "ft_multiplotER" took 1 seconds and an estimated 0 MB