[FieldTrip] Combined EEG and MEG source reconstruction
Mushfa Yousuf
mushfa.yousuf at googlemail.com
Fri Mar 8 16:12:52 CET 2013
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> 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
> http://mailman.science.ru.nl/mailman/listinfo/fieldtrip
>
>
>
> _______________________________________________
> fieldtrip mailing list
> fieldtrip at donders.ru.nl
> http://mailman.science.ru.nl/mailman/listinfo/fieldtrip
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20130308/2d03dd17/attachment.html>
More information about the fieldtrip
mailing list