[FieldTrip] Computing source-level ROI EEG connectivity

Schoffelen, J.M. (Jan Mathijs) jan.schoffelen at donders.ru.nl
Thu May 7 09:17:10 CEST 2020

Dear Michael,

It looks that you already got quite far! I agree that the DTF values do not make much sense.

a few pointers:

-I would use some regularisation for  the MNE reconstruction.
-I would check whether the individual trials after source reconstruction yield meaningful values in the ‘mom’-field. Particularly, are the different trials different? Do the time series look meaningful to begin with (both in terms of signal magnitude and ’shape’)? The reason I ask, is that I am not fully sure whether the ‘rawtrial’ option is behaving according to your expectations.
-Another way for getting the mom would be to use cfg.keepfilter in ft_sourceanalysis, and left-multiply each dipole location’s ’spatial filter’ (in sourceEEG.avg.filter) with a dataEEG.trial matrix.
-IMPORTANT: don’t use the ‘pow’ for the subsequent spectral analysis and connectivity estimation. You should work with the mom.
-In order to get a ‘parcellated’ representation of your source level data, I recommend to use ft_sourceparcellate, which should achieve the parcellation you performed by hand (with the binary masks).
-The parcellation operation can be an average of the mom across the dipoles that belong to a parcel, but an alternative might be the the 1st principal component. The advantage of the latter would be that the per-parcel signal will be reduced to a single time series (rather then 3 x/y/z orientation time courses), which is convenient (and required) for the interpretation of the subsequent computations.
-You may want to think a bit about the depth weighting, since the MNE favors superficial locations in terms of amplitude, and this will also percolate to the way the parcellated time series are computed.
-IMPORTANT: using ‘dtf’ without any additional specification of parameters is probably not what you want here. The reason is that under the hood a multivariate factorization of the spectral matrix is performed (at least: an attempt is made), which will miserably fail in your case. The reason for this failure (thinking a bit more about this, I am quite surprised that it worked to begin with) is the fact that you ask for a factorization of a size NparcelxNparcel(xfrequency bins) matrix, where the effective rank of the matrix is less than the number of channels in your EEG array. Even if Nparcel is less than the number of EEG electrodes, this factorization will be very poorly estimated and nonsensical. I’d recommend using ‘granger’ as a method, in combination with cfg.granger.sfmethod = ‘bivariate’. This results in a parcel-pairwise factorization (which is hopefully stable most of the time). I am not a fan of dtf/pdc anyway given that the measures are claimed to provide estimates unproblematic estimates that are interpretable without any caveats.
-ALSO IMPORTANT: if you go for a granger-like measure (granger/dtf/pdc), in the spectral analysis you should probably use a decent amount of zero padding (cfg.pad) and sufficient tapsmofrq in order to get spectral data that stands a chance to get a stable factorization.

So, no definitive answers unfortunately, but perhaps some food for thought that might get you closer to a good result.

Good luck and keep up the good work,


J.M.Schoffelen, MD PhD
Associate PI, VIDI-fellow - PI, language in interaction
Telephone: +31-24-3614793
Physical location: room 00.028
Donders Centre for Cognitive Neuroimaging, Nijmegen, The Netherlands

On 6 May 2020, at 22:34, Michael Glassen - Biomedical Engineer <MGlassen at kesslerfoundation.org<mailto:MGlassen at kesslerfoundation.org>> wrote:

Hi all,

I am having issues using fieldtrip for source reconstruction and connectivity analysis with preprocessed EEG data. I believe the issue is occurring when I try to project my sources to specific ROIs before getting the connectivity. I have outlined my pipeline below:

-Because I do not have subject MRIs, I am using a template headmodel and sourcemodel that come included with fieldtrip
               - headmodel = standard_bem.mat
               -sourcemodel = cortex_8196.surf.gii

-After importing EEG data, I project the EEG electrodes to the scalp using the code below:
               if strcmp(headmodel.type, 'bemcp')
             scalp_index = 3;
elseif strcmp(headmodel.type, 'dipoli')
             scalp_index = 1;
cfg = [];
cfg.method = 'project';
cfg.headshape = headmodel.bnd(scalp_index);

eegData.elec = ft_electroderealign(cfg, eegData.elec);

-I then calculate covariance on a trial basis and calculate the leadfield matrix:

               cfg = [];
cfg.covariance = 'yes';
cfg.covariancewindow = [-inf 0];
cfg.keeptrials = 'yes';
tlckEeg = ft_timelockanalysis(cfg, eegData);
cfg         = [];
cfg.elec    = tlckEeg.elec;   % sensor information
cfg.channel = tlckEeg.label;  % the used channels
cfg.grid    = sourcemodel;   % source points
cfg.headmodel = headmodel;   % volume conduction model
cfg.singleshell.batchsize = 5000; % speeds up the computation
leadfield   = ft_prepare_leadfield(cfg);

-Now with the leadfield calculated, I calculate sources on a trial basis here:

cfg               = [];
cfg.method        = 'mne';
cfg.sourcemodel   = leadfield;
cfg.headmodel     = headmodel;
cfg.mne.prewhiten = 'yes';
cfg.mne.lambda    = 0;
cfg.mne.scalesourcecov = 'yes';
cfg.rawtrial      = 'yes';
sourceEEG          = ft_sourceanalysis(cfg,tlckEeg);
               sourceEEG =

                              struct with fields:

                              time: [1×357 double]
                              inside: [8196×1 logical]
                              pos: [8196×3 double]
                               tri: [16384×3 double]
                              method: 'rawtrial'
                              trial: [1×400 struct]
                              df: 400
                              cfg: [1×1 struct]

               sourceEEG.trial =

                              1×400 struct array with fields:



ans =

                                             8196         357

-So at this point I have 400 trials, and for each of these trials I have a time course of source power at 8196 voxels and 357 time points for each of those voxels.
-I then import two atlases, one that includes broddman areas and one that includes left and right hemispheres. I interpolate these atlases onto my source model and get the following:
               broddman =
struct with fields:

                                      pos: [8196×3 double]
                                       tri: [16384×3 double]
                                        unit: 'm'
                                       tissue: [8196×1 double] (each of these is a value 1-70, corresponding to which atlas entry the voxel is located in)
                                       tissuelabel: {1×70 cell}(ROI labels for the atlas, includes broddman areas not separated by hemisphere)
                                      cfg: [1×1 struct]
hemis2 =
struct with fields:
pos: [8196×3 double]
                                      tri: [16384×3 double]
                                        unit: 'm'
                                        tissue: [8196×1 double] (each of these is a value 1-7, corresponding to which atlas entry the voxel is located in)

                                         tissuelabel: {1×7 cell} (ROI Labels for hemisphere atlas, contains 3 Left areas, 3 Right Areas and one inter-hemispheric)
                                      cfg: [1×1 struct]

Now I create a binary mask for each left and right broddman area. For each ROI in the broddman area atlas(1-70), I find all voxels that are located within. This gives me a list of voxels for each broddman area. I then take each of these lists and find which are in the left hemisphere and which are in the right to get left and right broddman areas.
After these steps I am left with 140 binary masks, (size 1x8196), which equal 1 if the voxel is contained in the specific ROI and 0 if it is not. Using these masks, I extract the time course power for each ROI by taking the mean power of every voxel in the roi at each time point:
                    roiData(i,i3,i2) = mean(sourceEEG.trial(i2).pow(Masks{i},i3));

where i = roi number,

       i2 = trial number,

       i3 = time point

In order to get this data into fieldtrip format, I put ROI data into an EEGLAB structure along with 6 emg channels that were recorded at the same time and use eeglab2fieldtrip function to have a fieldtrip structure containing the broddman area sources as channels. (contained in header)

After this I perform the final step of computing dtf connectivity between each of the source rois and the emg channels with the code below:

cfg           = [];
cfg.method    = 'mtmfft';
cfg.taper     = 'dpss';
cfg.output    = 'fourier';
cfg.tapsmofrq = 2;
freq          = ft_freqanalysis(cfg, header);
cfg           = [];
cfg.method    = 'dtf';
dtf           = ft_connectivityanalysis(cfg, freq);

The values I get in DTF are extremely low, 10^-20. I am not sure why but I believe it has to do with the way I am projecting source activations onto specific ROIs. Any help or advice on how to change my pipeline is appreciated.

Michael Glassen

  ­­   _______________________________________________
fieldtrip mailing list

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20200507/0f27c5b4/attachment-0001.htm>

More information about the fieldtrip mailing list