[FieldTrip] problem with EEG source reconstruction

Thomas Sauvigny tsauvigny at gmail.com
Sat Feb 22 14:32:02 CET 2014

Dear fellow fiedtrippers,

I would like to do a source reconstruction of EEG data using The FielTrip example MRI dataset via minimum-norm estimate. I stick close to the MNE tutorial (http://fieldtrip.fcdonders.nl/tutorial/minimumnormestimate) and use a few steps from the EEG headmodel tutorial as well (http://fieldtrip.fcdonders.nl/tutorial/headmodel_eeg). Unfortunately, I am struggling to combine the information from both tutorials in some crucial step prior to ft_sourceanalysis.
My first approach has been to basically run through the MNE tutorial using the Fieltrip MRI dataset. As I could not save both the original anatomy and the masked anatomy in a freesurfer compatible format while using 

cfg           = [];
cfg.output    = {'brain','skull','scalp'};
segmentedmri  = ft_volumesegment(cfg, mri);

I instead proceeded as mentioned below and aligned the electrodes to the scalp in a second step (according to the headmodel tutorial). But when I want to combine the headmodel (vol) and the electrodes (elec_aligned2) into the leadfield, that does not work:

using headmodel specified in the configuration
using electrodes specified in the configuration
Error using ft_prepare_vol_sens (line 534)
unsupported volume conductor model for EEG

Error in prepare_headmodel (line 94)
[vol, sens] = ft_prepare_vol_sens(vol, sens, 'channel', cfg.channel, 'order',

Error in ft_prepare_leadfield (line 124)
[vol, sens, cfg] = prepare_headmodel(cfg, data);

When I instead use the volume conduction model built in the headmodel tutorial (vol2) and the sourcespace from the MNE tutorial I get:

Error using  * 
Inner matrix dimensions must agree.

Error in minimumnormestimate (line 179)
      A       = P*A; % prewhitened leadfields

Error in ft_sourceanalysis (line 869)
        dip(i) = minimumnormestimate(grid, sens, vol, squeeze_avg, optarg{:},
        'noisecov', squeeze(Cy(i,:,:)));

In older posts in the mailing list there appeared to be problems with the scaling of the different components and I adapted the units, but that did not solve the problems.
I would be very thankful for any help.


My script:

clear all;
mri = ft_read_mri('Subject01.mri');
save mri mri;

cfg            = [];
cfg.resolution = 1;
cfg.dim        = [256 256 256];
mrirs          = ft_volumereslice(cfg, mri);
save mrirs mrirs;

% segmentation of the mri
% load mrirs;
cfg           = [];
cfg.coordsys  = 'ctf';
cfg.output    = {'skullstrip' 'brain'};
seg           = ft_volumesegment(cfg, mrirs);
save seg seg;


%%% Thus, you need to identify the anterior commissure, the posterior commissure, and an interhemispheric point (which defines the XZ-plane). You can navigate across the slices with the arrows. 
%%% By pressing the 'a', 'p', and 'z' at the right locations, the voxel
%%% coordinates of these landmarks are stored:
% load mrirs;
cfg        = [];
cfg.method = 'interactive';
mri_tal    = ft_volumerealign(cfg, mrirs);
save mri_tal mri_tal

clear all;
load mri_tal;
load seg;
% ensure that the skull-stripped anatomy is expressed in the same coordinate system as the anatomy
seg.transform = mri_tal.transform;
% save both the original anatomy, and the masked anatomy in a freesurfer compatible format
cfg             = [];
cfg.filename    = 'Subject01';
cfg.filetype = 'mgz';
cfg.parameter   = 'anatomy';
ft_volumewrite(cfg, mri_tal);
cfg.filename    = 'Subject01masked';
ft_volumewrite(cfg, seg);

%%%%% MAC Terminal Input %%%%%%%%%%%%%%%%
chsh -s /bin/tcsh
setenv FREESURFER_HOME /Applications/freesurfer
source $FREESURFER_HOME/SetUpFreeSurfer.csh

setenv SUBJECTS_DIR /Applications/freesurfer/subjects
setenv SUBJECT /Subject01
%%% .mgz files to mri folder
cd $SUBJECTS_DIR/Subject01/mri/
mri_convert -c -oc 0 0 0 Subject01masked.mgz orig.mgz
mri_convert -c -oc 0 0 0 Subject01.mgz orig-nomask.mgz

recon-all -talairach -subjid Subject01
recon-all -nuintensitycor -subjid Subject01
recon-all -normalization -subjid Subject01
cp T1.mgz brainmask.mgz
recon-all -gcareg -subjid Subject01
recon-all -canorm -subjid Subject01
recon-all -careg -subjid Subject01
recon-all -calabel -subjid Subject01
recon-all -normalization2 -subjid Subject01
recon-all -segmentation -subjid Subject01
recon-all -fill -subjid Subject01

%%% This ends the part of the Freesurfer pipeline concerned with volumetric processing. At this stage you should have a file filled.mgz containing the segmentation of the cortical white matter (cerebellum is not included!). 
%%% You can check how this looks using FieldTrip, by doing the following:
% % go to the Subject01/mri directory
mri = ft_read_mri('filled.mgz');
cfg = [];
cfg.interactive = 'yes';
figure;ft_sourceplot(cfg, mri);

%%% back in Terminal:
cd $SUBJECTS_DIR/Subject01/mri/
recon-all -tessellate -subjid Subject01
recon-all -smooth1 -subjid Subject01
recon-all -inflate1 -subjid Subject01
recon-all -qsphere -subjid Subject01
recon-all -fix -subjid Subject01
cp brain.mgz brain.finalsurfs.mgz
recon-all -finalsurfs -subjid Subject01
recon-all -smooth2 -subjid Subject01
recon-all -inflate2 -subjid Subject01
recon-all -sphere -subjid Subject01
recon-all -surfreg -subjid Subject01

%%% MNE download and installation %%%%%%%%%%%%%%
setenv MNE_ROOT /Applications/MNE-2.7.4-3378-MacOSX-x86_64
setenv MATLAB_ROOT /Applications/MATLAB_R2013b.app
source $MNE_ROOT/bin/mne_setup

%  http://martinos.org/mne/stable/manual/list.html
cd $MNE_ROOT/bin

setenv SUBJECTS_DIR /Applications/freesurfer/subjects
setenv SUBJECT /Subject01

./mne_setup_source_space --ico -6 --overwrite

%%%% back in Matlab: /Applications/freesurfer/subjects/Subject01/bem
bnd = ft_read_headshape('Subject01-oct-6-src.fif', 'format', 'mne_source');

%%% Current Folder : /Applications/freesurfer/subjects/Subject01/Subject01/mri
mri_nom = ft_read_mri('orig-nomask.mgz');

cfg = [];
cfg.method = 'interactive';
mri_nom_ctf = ft_volumerealign(cfg, mri_nom);

%%% fiducials chosen, fiducials l,r,n

mri_nom_ctf = ft_convert_units(mri_nom_ctf, 'cm');
T   = mri_nom_ctf.transform*inv(mri_nom_ctf.transformorig);
%%% Current Folder : /Applications/freesurfer/subjects/Subject01/Subject01/bem
bnd  = ft_read_headshape('Subject01-oct-6-src.fif', 'format', 'mne_source');
sourcespace = ft_convert_units(bnd, 'cm');
sourcespace = ft_transform_geometry(T, sourcespace);
save sourcespace sourcespace;
save T T; %we will need the transformation matrix in the next step

%%% Current Folder auf: /Applications/freesurfer/subjects/Subject01/Subject01/mri
mri_nom = ft_read_mri('orig-nomask.mgz');

cfg           = [];
cfg.coordsys  = 'spm'; 
cfg.output    = {'brain'};
seg           = ft_volumesegment(cfg, mri_nom);
seg           = ft_convert_units(seg,'cm');

cfg           = [];
cfg.method    = 'singleshell';
vol           = ft_prepare_headmodel(cfg,seg);
vol.bnd       = ft_transform_geometry(T, vol.bnd);
save vol vol;

% load vol                                       % volume conduction model
figure;hold on;
ft_plot_vol(vol, 'facecolor', 'none');alpha 0.5;
ft_plot_mesh(sourcespace, 'edgecolor', 'none'); camlight 

%%% so far, so good. >>> then I create the elec file:
clear all;
mri2 = ft_read_mri('Subject01.mri');
save mri2 mri2;

cfg           = [];
cfg.output    = {'brain','skull','scalp'};
segmentedmri2  = ft_volumesegment(cfg, mri2);
save segmentedmri2 segmentedmri2

cfg.numvertices = [3000 2000 1000];

save bnd2 bnd2


cfg        = [];
cfg.method ='dipoli';
vol2        = ft_prepare_headmodel(cfg, segmentedmri2);

save vol2 vol2


vol2 = ft_convert_units(vol2, 'cm');

ft_plot_mesh(vol2.bnd(1),'facecolor','none'); %scalp
ft_plot_mesh(vol2.bnd(2),'facecolor','none'); %skull
ft_plot_mesh(vol2.bnd(3),'facecolor','none'); %brain

ft_plot_mesh(vol2.bnd(1), 'facecolor',[0.2 0.2 0.2], 'facealpha', 0.3, 'edgecolor', [1 1 1], 'edgealpha', 0.05);
hold on;
hold on;
ft_plot_mesh(vol2.bnd(3),'edgecolor','none','facecolor',[0.4 0.6 0.4]);

elec2       = ft_read_sens('standard_1020.elc');   % may you need to define the path to the file


elec2 = ft_convert_units(elec2, 'cm');

% load volume conduction model
load vol2;                              
% head surface (scalp)
ft_plot_mesh(vol2.bnd(1), 'edgecolor','none','facealpha',0.8,'facecolor',[0.6 0.6 0.8]); 
hold on;
% electrodes
ft_plot_sens(elec2,'style', 'sk');    

mri2 = ft_convert_units(mri2, 'cm');

nas=warp_apply(transm,nas, 'homogenous');
lpa=warp_apply(transm,lpa, 'homogenous');
rpa=warp_apply(transm,rpa, 'homogenous');

fid.chanpos       = [nas; lpa; rpa];       % ctf-coordinates of fiducials
fid.label         = {'Nz','LPA','RPA'};    % same labels as in elec 
fid.unit          = 'cm';                  % same units as mri
% alignment
cfg               = [];
cfg.method        = 'fiducial';            
cfg.template      = fid;                   % see above
cfg.elec          = elec2;
cfg.fiducial      = {'Nz', 'LPA', 'RPA'};  % labels of fiducials in fid and in elec
elec_aligned2      = ft_electroderealign(cfg);

save elec_aligned2 elec_aligned2;
% elec_aligned2 = ft_convert_units(elec_aligned2, 'cm');
hold on;
ft_plot_mesh(vol2.bnd(1),'facealpha', 0.85, 'edgecolor', 'none', 'facecolor', [0.65 0.65 0.65]); %scalp

cfg           = [];
cfg.method    = 'interactive';
cfg.elec      = elec_aligned2;
cfg.headshape = vol2.bnd(1);
elec_aligned2  = ft_electroderealign(cfg);

%%%%% and here I am struggling: I can not get any of the above files to
%%%%% work together: vol or vol2, elec_aligned2, sourcespace etc.

vol = ft_convert_units(vol, 'cm');
elec_aligned2 = ft_convert_units(elec_aligned2, 'cm');
cfg = [];
cfg.elec      = elec_aligned2;                     % sensor positions
cfg.channel       = 'all';                   % the used channels
cfg.grid.pos = sourcespace.pnt;              % source points
cfg.grid.inside = 1:size(sourcespace.pnt,1); % all source points are inside of the brain
cfg.vol = vol;
cfg.normalize       = 'yes';% volume conduction model
leadfield = ft_prepare_leadfield(cfg);

save leadfield leadfield;      

clear all;
load leadfield;
load vol2;
load elec_aligned2;

cfg        = [];
cfg.method = 'mne';
cfg.elec = elec_aligned2;
cfg.grid   = leadfield;
cfg.vol    = vol;
cfg.mne.prewhiten = 'yes';
cfg.mne.lambda    = 3;
cfg.mne.scalesourcecov = 'yes';
cfg.mne.normalize = 'yes';
sourceCERT0  = ft_sourceanalysis(cfg, GA1);
sourceCERT3 = ft_sourceanalysis(cfg, GA2);
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20140222/1d1971f3/attachment-0001.html>

More information about the fieldtrip mailing list