[FieldTrip] Issue plotting Brainnetome aligned EEG source data

Schoffelen, J.M. (Jan Mathijs) janmathijs.schoffelen at donders.ru.nl
Thu Jul 1 20:17:55 CEST 2021


Hi Jack,

1. "It is a bit unclear to me what the relevance of the atlas is in all this."

Ultimately, I am looking to explore EEG connectivity between brainnetome ROIs. Specifically, in my full pipeline I am:
a) Creating Brainnetome atlas aligned template grid
b) Using the template grid for source reconstruction across all EEG trials (i.e. epochs)
c) Muliplying the spatial filter with the individual EEG epochs to get trial-level source data
d) Projecting the dipole moments to their max in each trial
e) Using ft_volumelookup to select and then average the dipole mom for each brainnetome ROI to create a 'virtual channel' timecourse
f) Running ft_freqanalysis and ft_connectivity analysis on the virtual channel data

Gotcha. I already assumed that overall this would be your goal. However, I would recommend a slightly different sequence of steps (involving some other functions):

a) interpolate the Brainnetome atlas to the template grid, using ft_sourceinterpolate. This requires the grid to be in the same coordinate system as the atlas, which is most likely the case if you use one of the standard_sourcemodel3d<some-number>mm.mat. This approach is slightly different from the recipe you used so far, which I realize comes from here: https://www.fieldtriptoolbox.org/example/create_single-subject_grids_in_individual_head_space_that_are_all_aligned_in_brain_atlas_based_mni_space/. This piece of documentation is a bit outdated in my opinion. It should be possible to do something like this:

atlas_grid = ft_sourceinterpolate(cfg, atlas, template_grid), using cfg.interpmethod = ’nearest’ (or ‘mode’)

b) good, don’t forget to use cfg.<method>.keepfilter = ‘yes’;

c),d),e) this triplet of steps is somewhat outdated. You should be able to wrap these steps into a single function call to ft_virtualchannel:

vcdata = ft_virtualchannel(cfg, data, source, atlas_grid); % please consult the documentation and the code of the function about how best to use this.

f) all good.


2. "The atlas should not be used a priori to define the ‘inside’ compartment. You can map the dipoles onto the anatomical parcels after you performed the source reconstruction."

Is this still the case if I am using a standard BEM model? I decided to create a template grid pre-aligned to the Brainnetome atlas to make it simple to compare subject source-level data and interpret the outcomes relative to specific ROIs in the brainnetome atlas. Is this a poor approach?

No, this is in principle a good approach. Note that you use a template for the volume conduction model, so I’d assume that the subject-specific grids are not different across subjects either.


I did first try the typical 'source reconstruction -> interpolation to atlas -> source parcellate' pipeline but I ran into a lot of problems trying to interpolate my source reconstruction onto the Brainnetome atlas. Further digging led me to the current approach as a possible solution.

As per what I wrote above I’d recommend to go the other way around: interpolate the atlas onto the lower resolution dipole grid. This makes life much easier.

3. "Following up on my previous message: I didn’t think of the fact that you are using EEG, with a BEM volume conduction model. In this case the inwardshift is probably not going to work."

You are right there, as far I can figure. I played a lot with various inwardshift, resolution, and other settings when I was preparing the template_grid (sourcemodel), and I also mucked around with ft_sourceplot projection methods. In some cases, certain options did look like they improved the surface plot but I think they were only hiding some of the problem. Inward shift only appeared to cover the surface when I had shifted many gridpoints 'outside' the brain volume; although, regardless, as soon as I align the template grid to the brainnetome atlas the problem returns.

OK, let’s forget about this then.

4. "Is there a particular reason why you’d want to visualize your results on the surface?"

Self-gratification, perhaps. Really, I was hoping to plot the outcomes on a surface because it is a nice and simple way to display results, both for checking that things 'look sensible' and for potential publications where we might use this atlas. Is
it the data type that makes this tricky? I am wondering if my pipeline is ok, but perhaps I need to give up on the surface plot, or create one specifically to represent data at ROIs of the Brainnetome atlas.

In case it is helpful for you, or for other interested readers, my full code for steps a-e is below. Thanks again for your thoughts!


I hope the above pointers are helpful in getting you on the road again.

Best wishes,
Jan-Mathijs



Best,
Jack



%% Fieldtrip: eLORETA ROI timecourses
% https://urldefense.com/v3/__https://www.fieldtriptoolbox.org/example/create_single-subject_grids_in_individual_head_space_that_are_all_aligned_in_brain_atlas_based_mni_space/__;!!HJOPV4FYYWzcc1jazlU!q8K-aqjtprbNX-bLhVdFRIr1lKXtT8VKw6LQhGjxl13kWaxe7G9gYp1JvXT18OEIpcZlnTLVcH9McUU$
% https://urldefense.com/v3/__https://www.fieldtriptoolbox.org/tutorial/salzburg/*reconstruct-the-time-course-of-activity-at-a-particular-brain-location__;Iw!!HJOPV4FYYWzcc1jazlU!q8K-aqjtprbNX-bLhVdFRIr1lKXtT8VKw6LQhGjxl13kWaxe7G9gYp1JvXT18OEIpcZlnTLVl8J_GVQ$

   ft_defaults
   [ftver, ftpath] = ft_version;

   % Load the Headmodel
   load(fullfile(ftpath, 'template/headmodel/standard_bem.mat')); % Colin27

   % Prepare the template grid 'standard sourcemodel'
   cfg = [];
   cfg.xgrid  = 'auto';
   cfg.ygrid  = 'auto';
   cfg.zgrid  = 'auto';
   cfg.unit   = 'mm';
   cfg.tight  = 'yes';
   cfg.resolution  = 7; % What does this result in regard to the final resolution for reporting?
   cfg.headmodel   = vol;
   template_grid   = ft_prepare_sourcemodel(cfg);

   % Check alignment of the template grid and headmodel
   figure;
   hold on
   ft_plot_mesh(template_grid.pos(template_grid.inside,:));
   ft_plot_headmodel(vol,  'facecolor', 'cortex', 'edgecolor', 'none');
   ft_plot_axes(vol);
   alpha 0.5
   camlight

   % Read in Brainnetome Atlas; ensure its units are consistent with template
   atlas = ft_read_atlas(fullfile(ftpath, 'template/atlas/brainnetome/BNA_MPM_thr25_1.25mm.nii'));

   % Ensure units are consistent between atlas and sourcemodel
   atlas = ft_convert_units(atlas,'mm');
   template_grid = ft_convert_units(template_grid,'mm');

   % Check the coordsys field exists in sourcemodel
   if ~isfield(template_grid, 'coordsys');
       template_grid.coordsys = 'mni';
   end

   % Create binary mask for all grid points inside the Brainnetome atlas locations
   cfg = [];
   cfg.atlas      = atlas;
   cfg.roi        = atlas.tissuelabel;  % here you can also specify a single label, i.e. single ROI
   cfg.inputcoord = 'mni';
   mask           = ft_volumelookup(cfg, template_grid);
   template_grid.inside = false(template_grid.dim);
   template_grid.inside(mask == 1) = true;

   % Check grid points inside the Brainnetome atlas
   figure;
   ft_plot_mesh(template_grid.pos(template_grid.inside,:),'vertexcolor','black','vertexalpha',0.5);

   clearvars -except atlas template_grid vol

   %% Prep EEG source analysis

   setenv('PATH', 'C:\Program Files\OpenMEEG\bin');
   setenv('LD_LIBRARY_PATH', 'C:\Program Files\OpenMEEG\lib');

   % Load preprocessed EEG data
   load('C:\Users\jf752\Desktop\Processing\source_analysis\test_data');

   % Realign electrodes
   cfg = [];
   cfg.method = 'interactive';
   cfg.headshape = vol.bnd(1);
   cfg.elec = cln_data.elec;
   elec_aligned = ft_electroderealign(cfg);

   cln_data.elec = elec_aligned;

   % Convert time for each trial so that they are equivalent (for rawtimelock)
   for i = 1:length(cln_data.time);
       cln_data.time{1,i} = 0:0.002:1.998;
   end

   clear i elec_aligned cfg

   % Compute sensor-level time-locked data
   cfg = [];
   cfg.trials = 'all';
   cfg.keeptrials = 'yes';
   cfg.covariance = 'yes';

   timlock = ft_timelockanalysis(cfg,cln_data);

   % Prepare the leadfield
   cfg = [];
   cfg.elec      = cln_data.elec;
   cfg.channel   = cln_data.label;
   cfg.headmodel = vol; % volume conduction model
   cfg.sourcemodel.pos    = template_grid.pos; % number of voxels

   leadfield = ft_prepare_leadfield(cfg);

   clear cfg

   %% Run source localisation and calculate single-trial ROI timecourses

   % Compute eLORETA time-course based on the total data
   cfg = [];
   cfg.method       = 'eloreta';
   cfg.elec         = timlock.elec;
   cfg.sourcemodel  = leadfield;
   cfg.headmodel    = vol;
   cfg.projectnoise = 'yes';
   cfg.keepmom      = 'yes';
   cfg.keepfilter   = 'yes';

   source = ft_sourceanalysis(cfg, timlock);

   % Calculate 3D trial moms based on Michael Glassen's FTMailing Request
   % https://mailman.science.ru.nl/pipermail/fieldtrip/2020-May/039989.html
   for f = 1:length(source.avg.filter);
       for t = 1:length(cln_data.trial);
           if isempty(source.avg.filter{f});
               x = [];
           else
               x = source.avg.filter{f}*cln_data.trial{t};
           end
       trial_source.trial(t).mom{f} = x;
       end
   end

   % Mimic Source.avg config to project trial dipole moments in direction of maximal power
   tmp = source;

   for t = 1:length(cln_data.trial);
       tmp.avg.mom = trial_source.trial(t).mom;

       cfg = [];
       cfg.projectmom = 'yes';
       tmp = ft_sourcedescriptives(cfg, tmp);

       trial_source.trial(t).mom = tmp.avg.mom;
       trial_source.trial(t).pow = tmp.avg.pow;
       trial_source.trial(t).ori = tmp.avg.ori;
       trial_source.trial(t).label = tmp.avg.label;
       trial_source.trial(t).filterdimord = tmp.avg.filterdimord;
   end

   clear f t x tmp ans cfg

   % Calculate the average signal in each ROI
   for t = 1:length(cln_data.trial);
       for r = 1:length(atlas.tissuelabel);
           cfg = [];
           cfg.atlas      = atlas;
           cfg.roi        = atlas.tissuelabel{1,r};  % here you can also specify a single label, i.e. single ROI
           cfg.inputcoord = 'mni';
           mask           = ft_volumelookup(cfg, template_grid);

           mask_idx = find(mask);

           for m = 1:length(mask_idx);
               if isempty(trial_source.trial(t).mom{mask_idx(m)});
                  tmp(m,:) = nan(1,length(cln_data.trial{1}));
               else
                  tmp(m,:) = trial_source.trial(t).mom{mask_idx(m)};
               end
           end

           roi(r,:) = nanmean(tmp);
           clear tmp mask mask_idx ans m
       end
       roi_trial_dat{t} = roi;
       clear roi r
   end

_______________________________________________
fieldtrip mailing list
https://mailman.science.ru.nl/mailman/listinfo/fieldtrip
https://urldefense.com/v3/__https://doi.org/10.1371/journal.pcbi.1002202__;!!HJOPV4FYYWzcc1jazlU!q8K-aqjtprbNX-bLhVdFRIr1lKXtT8VKw6LQhGjxl13kWaxe7G9gYp1JvXT18OEIpcZlnTLVNDcHIOs$

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20210701/26cd5080/attachment.htm>


More information about the fieldtrip mailing list