[FieldTrip] Combined EEG and MEG source reconstruction

"Jörn M. Horschig" jm.horschig at donders.ru.nl
Mon Mar 11 10:41:33 CET 2013


whoopsy, just realized that you use singlesphere, not singleshell. I was 
talking about singleshell ;)

On 3/11/2013 9:47 AM, "Jörn M. Horschig" wrote:
> Dear Mushfa,
>
> I never tried what you did, however from what I know of the matter and 
> from what Robert wrote, I can see that you made some errors in your 
> script that I can point to. I can, however, neither guarantee that 
> this solves all or any of your problems nor that it's all correct what 
> you are doing. So, I can only tell you some  minor things:
> 1) When defining trial, setting cfg.trialdef.eventtype to '?' should 
> actually only be done if you do not know the eventtype. Also note that 
> it might be better to use all your trials rather than a subset - 
> better signal-to-noise that way ;)
>
> 2) No reason to segment the brain twice - this operation should be 
> done only once and the segmented brain you base the two 
> leadfield-matrices on *must* be exactly the same. So you can save time 
> by removing the computation of vol2.
>
> 3) Specifying only one frequency might not be a good idea, because of 
> imperfect sampling.
>
> 4) 'singlesphere' is created specifically for MEG data. BEM models are 
> good for EEG, not sure how good they are for MEG, but the dipoli 
> implementation only works under unix (Linux/Mac) and openMEEG has to 
> be downloaded manually first to make it work.
>
> 5) As Robert said, units are not always properly controlled for using 
> FieldTrip. You can use ft_convert_units to make them match.
>
> 6) Please check again how 'grid' is defined, e.g. the .dim field 
> shouldn't be concatenated.
>
> 7) You do not only need to concatenate the grids, but also the data, 
> hope you did that ;) But I would really rather concatenated the two 
> freq structures you used to compute the sourcemodel/leadfield than to 
> compute it newly. If you concatenate you can be sure that the channel 
> order is persistent.
>
>
> If you keep on having problems, I would suggest to use breakpoints or 
> /dbstop if error/ and then step through the function to see where the 
> error comes from, i.e. why there are NaNs or Infs in (lf*Cf*lf), and 
> start off with looking at if there are NaNs or Infs in lf or Cf 
> already. That's all I can help you with, cause that's were my 
> knowledge ends ;) Maybe someone else is wiser than I am. Good luck!
>
> Best,
> Jörn
>
> On 3/8/2013 4:12 PM, Mushfa Yousuf wrote:
>>
>> Hello;
>>
>> I have few problems while calculating combine "MEG & EEG" source 
>> reconstruction and I need your help to troubleshoot the problems.
>>
>> I have a data from a Parkinson's disease patient who underwent DBS 
>> and we recorded 306 MEG channels and 64 EEG channels simultaneously, 
>> first I have calculated the leadfield individually for both EEG and MEG.
>>
>> I have defined all the calculations step by step so that I can 
>> describe my problems very precisely.
>>
>> *_CALCULATE THE LEADFIELD FOR MEG CHANNELS:-_*
>>
>> // I have defined the Pre-Stimulus data i.e. Stim-off data in the 
>> original code.
>>
>> _Reading raw data set for MEG channels_
>>
>> __
>>
>> cfg = [];
>>
>> cfg.dataset  = 'ap_rest_stimon_bi_2_5ma.fif';
>>
>> cfg.trialdef.eventtype      = '?';
>>
>> cfg.trialdef.triallength = 1;
>>
>> cfg.trialdef.poststim       = 88;
>>
>> cfg = ft_definetrial(cfg);
>>
>> cfg.trl = cfg.trl(1:88,:);   // Read 88 trials
>>
>> cfg.channel   = {'MEG'}; *// In order to read raw data only for MEG 
>> channels*
>>
>> dataFIC12 = ft_preprocessing(cfg);
>>
>> After that I have defined the _time of window of interest_ for the 
>> MEG raw data set
>>
>> cfg.toilim = [0 0.5];
>>
>> dataPost = ft_redefinetrial(cfg, dataFIC12);
>>
>> and calculate the _Power Spectral density matrix_ for the MEG 
>> channels data set
>>
>> cfg = [];
>>
>> cfg.method    = 'mtmfft';
>>
>> cfg.output    = 'powandcsd';
>>
>> cfg.tapsmofrq = 2;
>>
>> cfg.foilim    = [129 129]; // Frequency of interest is 129 Hz the 
>> stimulation frequency
>>
>> cfg.pad = 1;
>>
>> freqPost = ft_freqanalysis(cfg, dataPost);
>>
>> Now next step for me is to choose the method for the Head Model. I 
>> have noticed that Singlesphere is the common method which can be used 
>> for both EEG and MEG. So I have selected the singlesphere method for 
>> the calculation of volume conduction model.
>>
>> mri = ft_read_mri('ms1879865-0004-00001-000192-01.img'); // READ MRI
>>
>> cfg = [];
>>
>> cfg.write      = 'no';
>>
>> cfg.coordsys   = 'ctf';
>>
>> segmentedmri1 = ft_volumesegment(cfg, mri); // Calculate the 
>> segmented mri
>>
>> cfg = [];
>>
>> cfg.method = 'singlesphere'; // *METHOD*
>>
>> vol1 = ft_prepare_headmodel(cfg, segmentedmri1);
>>
>> Now next step is to calculate the LEADFIELD
>>
>> cfg                 = [];
>>
>> cfg.grad            = freqPost.grad;
>>
>> cfg.vol             = vol1;
>>
>> cfg.reducerank      = 2;    // MEG
>>
>> cfg.channel         = {'MEG'};
>>
>> cfg.grid.resolution = 1; // I choose a resolution 1
>>
>> [grid1] = ft_prepare_leadfield(cfg);
>>
>> *disp(grid1)*
>>
>> *        xgrid: [-7 -6 -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7]*
>>
>> *        ygrid: [-7 -6 -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7 8]*
>>
>> *        zgrid: [-4 -3 -2 -1 0 1 2 3 4 5 6 7 8 9 10 11]*
>>
>> *          dim: [15 16 16]*
>>
>> *          pos: [3840x3 double]*
>>
>> *         unit: 'cm'*
>>
>> *       inside: [1x2037 double]*
>>
>> *      outside: [1x1803 double]*
>>
>> *          cfg: [1x1 struct]*
>>
>> *    leadfield: {1x3840 cell}*
>>
>> **
>>
>> **
>>
>> ØIt has to be noticed that Inside Dipole is 2037.
>>
>> **
>>
>> *_CALCULATE THE LEADFIELD FOR EEG CHANNELS:-_*
>>
>> _Reading raw data set for EEG channels_
>>
>> __
>>
>> cfg = [];
>>
>> cfg.dataset  = 'ap_rest_stimon_bi_2_5ma.fif';
>>
>> cfg.trialdef.eventtype      = '?';
>>
>> cfg.trialdef.triallength = 1;
>>
>> cfg.trialdef.poststim       = 88;
>>
>> cfg = ft_definetrial(cfg);
>>
>> cfg.trl = cfg.trl(1:88,:);   // Read 88 trials
>>
>> cfg.channel   = {'EEG'}; *// In order to read raw data only for EEG 
>> channels*
>>
>> dataFIC12 = ft_preprocessing(cfg);
>>
>> After that I have defined the _time of window of interest_ for the 
>> EEG raw data set
>>
>> cfg.toilim = [0 0.5];
>>
>> dataPost = ft_redefinetrial(cfg, dataFIC12);
>>
>> and calculate the _Power Spectral density matrix_ for the EEG 
>> channels data set
>>
>> cfg = [];
>>
>> cfg.method    = 'mtmfft';
>>
>> cfg.output    = 'powandcsd';
>>
>> cfg.tapsmofrq = 2;
>>
>> cfg.foilim    = [129 129]; // Frequency of interest is 129 Hz
>>
>> cfg.pad = 1;
>>
>> freqPost = ft_freqanalysis(cfg, dataPost);
>>
>> _HEAD MODEL FOR EEG_
>>
>> mri = ft_read_mri('ms1879865-0004-00001-000192-01.img');
>>
>> cfg = [];
>>
>> cfg.write      = 'no';
>>
>> cfg.coordsys   = 'ctf';
>>
>> segmentedmri2 = ft_volumesegment(cfg, mri);
>>
>> cfg = [];
>>
>> cfg.method = 'singlesphere';  // METHOD
>>
>> vol2 = ft_prepare_headmodel(cfg, segmentedmri2);
>>
>> _LEADFIELD:___
>>
>> cfg            = [];
>>
>> cfg.elec       = freqPost.elec;
>>
>> cfg.vol        = vol2;
>>
>> cfg.reducerank = 3; % 3 for EEG
>>
>> cfg.channel    = {'EEG'};
>>
>> cfg.grid.resolution = 9.54; // Increase the resolution so that the 
>> inside dipoles for
>>
>> EEG should be close to the inside dipoles of MEG
>>
>> [grid2] = ft_prepare_leadfield(cfg);
>>
>> *disp(grid2)*
>>
>> *        xgrid: [1x16 double]*
>>
>> *        ygrid: [1x16 double]*
>>
>> *        zgrid: [1x13 double]*
>>
>> *          dim: [16 16 13]*
>>
>> *          pos: [3328x3 double]*
>>
>> *         unit: 'mm'*
>>
>> *       inside: [1x2038 double]*
>>
>> *      outside: [1x1290 double]*
>>
>> *          cfg: [1x1 struct]*
>>
>> *    leadfield: {1x3328 cell}*
>>
>> ØIt has to be noticed that Inside Dipole is 2038.
>>
>> *_Concatenate of MEG and EEG LEADFIELD:-_*
>>
>> Now I have the individual LEADFIELD of both MEG and EEG grid1 and 
>> grid2 respectively.
>>
>> *disp(grid1)*
>>
>> *        xgrid: [-7 -6 -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7]*
>>
>> *        ygrid: [-7 -6 -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7 8]*
>>
>> *        zgrid: [-4 -3 -2 -1 0 1 2 3 4 5 6 7 8 9 10 11]*
>>
>> *          dim: [15 16 16]*
>>
>> *          pos: [3840x3 double]*
>>
>> *         unit: 'cm'*
>>
>> *       inside: [1x2037 double]*
>>
>> *      outside: [1x1803 double]*
>>
>> *          cfg: [1x1 struct]*
>>
>> *    leadfield: {1x3840 cell}*
>>
>> **
>>
>> *disp(grid2)*
>>
>> *        xgrid: [1x16 double]*
>>
>> *        ygrid: [1x16 double]*
>>
>> *        zgrid: [1x13 double]*
>>
>> *          dim: [16 16 13]*
>>
>> *          pos: [3328x3 double]*
>>
>> *         unit: 'mm'*
>>
>> *       inside: [1x2038 double]*
>>
>> *      outside: [1x1290 double]*
>>
>> *          cfg: [1x1 struct]*
>>
>> *    leadfield: {1x3328 cell}*
>>
>> **
>>
>> *Few things to be noticed here *
>>
>> **
>>
>> Ø*unit is different for both MEG 'cm' and EEG 'mm'*
>>
>> Ø*I have used the same method for calculating the Volume for both MEG 
>> & EEG => "Singlesphere" *
>>
>> Ø*Concatenating the "pos" vector of EEG and MEG and the other vectors*
>>
>> **
>>
>> *Grid.pos = cat(1,grid1.pos,grid2.pos);*
>>
>> *Grid.xgrid = [grid1.xgrid grid2.xgrid];*
>>
>> *Grid.ygrid = [grid1.ygrid grid2.ygrid];*
>>
>> *Grid.zgrid = [grid1.zgrid grid2.zgrid];*
>>
>> *Grid.inside = [grid1.inside grid2.inside];*
>>
>> *Grid.outside = [grid1.outside grid2.outside];*
>>
>> *Grid.leadfield = [grid1.leadfield grid2.leadfield];*
>>
>> *Grid.dim = [grid1.dim grid2.dim];*
>>
>> **
>>
>> *_Combined Source Reconstruction:-_*
>>
>> cfg = [];
>>
>> cfg.dataset  = 'ap_rest_stimon_bi_2_5ma.fif';
>>
>> cfg.trialdef.eventtype      = '?';
>>
>> cfg.trialdef.triallength = 1;
>>
>> cfg.trialdef.poststim       = 88;
>>
>> cfg = ft_definetrial(cfg);
>>
>> cfg.trl = cfg.trl(1:88,:);   // Read 88 trials
>>
>> cfg.channel   = {'MEG','EEG'}; *// In order to read raw data for MEG 
>> & EEG channels*
>>
>> dataFIC12 = ft_preprocessing(cfg);
>>
>> After that as previously defined, I have calculated the Spectral 
>> density.
>>
>> Now I have to calculate the inverse solution for both EEG & MEG data 
>> trials.
>>
>> _Source Analysis: without contrasting condition:-_
>>
>> __
>>
>> cfg = [];
>>
>> cfg.method = 'dics';
>>
>> cfg.frequency = 129;
>>
>> cfg.grid = Grid;  // Concatenated grid
>>
>> cfg.vol = vol1;  // As vol1 = vol2 : both singlesphere
>>
>> cfg.powmethod = 'trace';
>>
>> cfg.dics.projectnoise = 'yes';
>>
>> cfg.dics.lambda = 0.0001;
>>
>> sourcePost_nocon = ft_sourceanalysis(cfg, freqPost); // *freqpost 
>> contains both MEG *
>>
>> *and EEG row data*
>>
>> ****
>>
>> I have received the following error:
>>
>> scanning grid 2037/4075
>>
>> ??? Error using ==> svd
>>
>> Input to SVD must not contain NaN or Inf.
>>
>> Error in ==> beamformer_dics>pinv at 568
>>
>> [U,S,V] = svd(A,0);
>>
>> Error in ==> beamformer_dics at 314
>>
>> filt = pinv(lf' * invCf * lf) * lf' * invCf;              % Gross eqn. 3,
>>
>> use PINV/SVD to cover rank deficient leadfield
>>
>> Error in ==> ft_sourceanalysis at 595
>>
>> dip(i) = beamformer_dics(grid, sens, vol, [], squeeze(Cf(i,:,:)),
>>
>> optarg{:});
>>
>> Error in ==> code at 227
>>
>> sourcePost_nocon = ft_sourceanalysis(cfg, freqPost);
>>
>> ØIf you notice that till 2037 grid it represents the MEG grid and 
>> after that EEG grid is concatenated. i.e 2037 +2038 = 4075 (which 
>> also represent inside dipole).
>>
>> ØSecondly I try to use other method for the head modeling i.e 
>> openmeeg, dipoli and bemcp method which leads me to different errors 
>> which I was unable to trouble shoot.
>>
>> ·Openmeeg = OpenMEEG binaries are not correctly installed
>>
>> ·Dipoli =vol.mat file was missing after calculating the headmodel
>>
>> ·Bemcp =vol.mat file contains NaN fiels which leads to wrong
>>
>> leadfield computation.
>>
>> *_I apologize for the long email. I wanted to show the steps which I 
>> used so that it is easier to debug the problem._*
>>
>> *__*
>>
>> Thank you;
>>
>> Mushfa Yousuf
>>
>>
>>
>> On Fri, Feb 8, 2013 at 2:36 PM, Robert Oostenveld 
>> <r.oostenveld at donders.ru.nl <mailto:r.oostenveld at donders.ru.nl>> wrote:
>>
>>     Dear Mushfa
>>
>>     In principle it can be made to work, but you'll have to do a
>>     little bit of  work outside fieldtrip. You would start with
>>     1) ft_prepare_headmodel and ft_prepare_leadfield for the MEG part
>>     2) ft_prepare_headmodel and ft_prepare_leadfield for the EEG part
>>     Subsequently you would (manually) concatenate the leadfields of
>>     the MEG and EEG. Here you have to make sure that the channel
>>     order is consistent with that in your combined data.
>>     The resulting structure with combined leadfields can be passed
>>     into ft_sourceanalysis as cfg.grid.
>>
>>     However, I foresee a problem with the scaling. At this moment
>>     FieldTrip is not (yet) very precise in maintaining the correct
>>     units throughout. Consequently you might get different units for
>>     the EEG and MEG contribution to each dipole, causing weird
>>     cancelation effects. You'll just have to check carefully that the
>>     leadfield units and the data units consistently scaled between
>>     megmag, megplanar and eeg.
>>
>>     This brings me to another known problem for the combined
>>     vectorview magnetometer and planar gradiometer channels: at this
>>     moment the forward computation computes both in T (or fT),
>>     whereas the data is represented as T for the magnetometers and
>>     T/m for the gradiometers. So the gradiometer leadfields are off
>>     by a factor "1/baseline", which is approximately 1/70 if I recall
>>     correctly. Many neuromag users use only the planar channels, and
>>     therefore don't notice. We are working on fixing the neuromag
>>     channel/gradiometer units. That is not as simple as it appears,
>>     because conversions of the geometrical dimensions from m to cm or
>>     mm also affect the gradiometer baseline. So converting the data
>>     from m to mm should result in the planar channel data getting
>>     1000 times smaller, whereas the magnetometer channels shoudl stay
>>     the same. See http://bugzilla.fcdonders.nl/show_bug.cgi?id=963
>>     and http://bugzilla.fcdonders.nl/show_bug.cgi?id=963 if you want
>>     to follow this.
>>
>>     best regards
>>     Robert
>>
>>
>>
>>     On 28 Jan 2013, at 20:20, Mushfa Yousuf wrote:
>>
>>>     Hello Fieldtrippers,
>>>
>>>     I am currently doing source analysis using beamformer method in fieldtrip.
>>>
>>>     The data i work on is a simultaneous measurement of MEG (306) with EEG
>>>     (128) on PD patients.
>>>
>>>     I was successful in calculating the sources for both the modalities alone.
>>>
>>>     I tried according the tutorial "Combined EEG and MEG source reconstruction" and was successful till the estimation
>>>
>>>
>>>     of the lead field. But, then when i want to use the created vol for the
>>>     source analysis i have two structures one for EEG and the second one for
>>>     MEG. Now, i can only define only one to estimate the sources. How to use
>>>
>>>
>>>     this for both the modalities? Any help would be appreciated.
>>>
>>>     With regards
>>>     Mushfa
>>>     _______________________________________________
>>>     fieldtrip mailing list
>>>     fieldtrip at donders.ru.nl <mailto:fieldtrip at donders.ru.nl>
>>>     http://mailman.science.ru.nl/mailman/listinfo/fieldtrip
>>
>>
>>     _______________________________________________
>>     fieldtrip mailing list
>>     fieldtrip at donders.ru.nl <mailto:fieldtrip at donders.ru.nl>
>>     http://mailman.science.ru.nl/mailman/listinfo/fieldtrip
>>
>>
>>
>>
>> _______________________________________________
>> fieldtrip mailing list
>> fieldtrip at donders.ru.nl
>> http://mailman.science.ru.nl/mailman/listinfo/fieldtrip
>
>
> -- 
> Jörn M. Horschig
> PhD Student
> Donders Institute for Brain, Cognition and Behaviour
> Centre for Cognitive Neuroimaging
> Radboud University Nijmegen
> Neuronal Oscillations Group
> FieldTrip Development Team
>
> P.O. Box 9101
> NL-6500 HB Nijmegen
> The Netherlands
>
> Contact:
> E-Mail:jm.horschig at donders.ru.nl
> Tel:    +31-(0)24-36-68493
> Web:http://www.ru.nl/donders
>
> Visiting address:
> Trigon, room 2.30
> Kapittelweg 29
> NL-6525 EN Nijmegen
> The Netherlands
>
>
> _______________________________________________
> fieldtrip mailing list
> fieldtrip at donders.ru.nl
> http://mailman.science.ru.nl/mailman/listinfo/fieldtrip


-- 
Jörn M. Horschig
PhD Student
Donders Institute for Brain, Cognition and Behaviour
Centre for Cognitive Neuroimaging
Radboud University Nijmegen
Neuronal Oscillations Group
FieldTrip Development Team

P.O. Box 9101
NL-6500 HB Nijmegen
The Netherlands

Contact:
E-Mail: jm.horschig at donders.ru.nl
Tel:    +31-(0)24-36-68493
Web: http://www.ru.nl/donders

Visiting address:
Trigon, room 2.30
Kapittelweg 29
NL-6525 EN Nijmegen
The Netherlands

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20130311/78fd7b38/attachment-0001.html>


More information about the fieldtrip mailing list