[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