[FieldTrip] MEG connectivity analysis - artifacts in the connectivity matrix

[Kinga Korek]

Hello,

I've been trying to put together all the pieces of code I needed for calculating whole brain functional connectivity from MEG recordings. The code is running and gives some results, but the connectivity matrix it produces seems to have some artefacts. Could you advise what should I change to fix it?
I'm using data from the CamCAN project, and for now, I'm trying with a single subject.
Here is the script I'm using:

ft_defaults

meg_file = 'sub-CC110033_ses-rest_task-rest_proc-sss.fif';
mri_file = 'sub-CC110033_T1w.nii';

% from: https://natmeg.se/MEEG_course2018/preprocess_MRI_data2.html

%% HEAD MODEL %%

%% Reading files and determining coordinate system

mri = ft_read_mri(mri_file);

%% Reading Anatomical Landmark Coordinates from .jason file (MEG data folder)

fileName = 'sub-CC110033_ses-rest_task-rest_proc-sss_coordsystem.json'; % filename in JSON extension
fid = fopen(fileName); % Opening the file
raw = fread(fid,inf); % Reading the contents
str = char(raw'); % Transformation
fclose(fid); % Closing the file
coord_data = jsondecode(str); % Using the jsondecode function to parse JSON from string

NAS = coord_data.AnatomicalLandmarkCoordinates.NAS;
LPA = coord_data.AnatomicalLandmarkCoordinates.LPA;
RPA = coord_data.AnatomicalLandmarkCoordinates.RPA;

%% Co-register MRI to Neuromag based on fiducials

cfg = [];
cfg.method = 'fiducial';
cfg.fiducial.nas = NAS;
cfg.fiducial.lpa = LPA;
cfg.fiducial.rpa = RPA;
cfg.coordsys = 'acpc'; % only with acpc coordinate system co-registration looks alright

mri_realigned_1 = ft_volumerealign(cfg, mri);

headshape = ft_read_headshape(meg_file); % Load head shape
grad = ft_read_sens(meg_file,'senstype','meg'); % Load MEG sensors

%% Automatic procedure called “iterative closest point” (icp) to fit fiducials and head points
cfg = [];
cfg.method              = 'headshape';
cfg.headshape.headshape = headshape;
cfg.headshape.icp       = 'yes';
cfg.coordsys            = 'acpc';
cfg.headshape.interactive = 'no';

mri_realigned_2 = ft_volumerealign(cfg, mri_realigned_1);
% sensors do not fit the head shape at this point!

%% Checking co-registration

cfg = [];
cfg.method              = 'headshape';
cfg.headshape.headshape = headshape;
cfg.coordsys            = 'acpc';
cfg.headshape.icp       = 'no';
cfg.headshape.interactive = 'no';

mri_realigned_3 = ft_volumerealign(cfg, mri_realigned_2);
% sensors fit the head shape now

%% Reslice MRI data and convert units

cfg = [];
cfg.resolution = 1;
mri_resliced_mm = ft_volumereslice(cfg, mri_realigned_3);
mri_resliced_cm = ft_convert_units(mri_resliced_mm, 'cm');

%% Segment data into brain, skull, and scalp

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

%% Correct the segmentation

binary_brain = mri_segmented.brain;
binary_skull = mri_segmented.skull | binary_brain;
binary_scalp = mri_segmented.scalp | binary_brain | binary_skull;
binary_scalp = mri_segmented.scalp + binary_skull;

% use boolean logic together with IMERODE
binary_skull = binary_skull & imerode(binary_scalp, strel_bol(2)); % fully contained inside eroded scalp
binary_brain = binary_brain & imerode(binary_skull, strel_bol(2)); % fully contained inside eroded skull

%% Copy the segmented volume into a new structure

mri_segmented_2 = mri_segmented;  % Copy stucture
mri_segmented_2.anatomy = mri_resliced_cm.anatomy; % Copy anatomical image

% insert the updated binary volumes, taking out the center part for skull and scalp
mri_segmented_2.brain    = binary_brain;
mri_segmented_2.skull    = binary_skull & ~binary_brain;
mri_segmented_2.scalp    = binary_scalp & ~binary_brain & ~binary_skull;

%% Construct mesh for brain segment

cfg = [];
cfg.method = 'projectmesh';
cfg.tissue = 'brain';
cfg.numvertices = 3000;

mesh_brain = ft_prepare_mesh(cfg, mri_segmented_2);

%% Construct head model for MEG

cfg = [];
cfg.method = 'singleshell';

headmodel = ft_prepare_headmodel(cfg, mesh_brain);

% from: https://natmeg.se/MEEG_course2018/preprocess_MRI_data2.html
% END of NatMEG tutorial


%% SOURCE MODEL %%

%% Load source model template
ftpath = 'fieldtrip-20220707\template\sourcemodel\standard_sourcemodel3d10mm.mat';
load(ftpath);
template_grid = sourcemodel;
clear sourcemodel;

%% Prepare source model based on MNI template

cfg = [];
% cfg.warpmni = 'yes'; % Warning: please specify cfg.method='basedonmni' instead of cfg.warpmni='yes'
cfg.method = 'basedonmni';
cfg.nonlinear = 'yes';
cfg.unit = 'cm';
cfg.resolution = 1;
cfg.template = template_grid;
cfg.headshape = headshape;
cfg.headmodel = headmodel;
cfg.mri = mri_resliced_cm;

sourcemodel = ft_prepare_sourcemodel(cfg);

%% ########################################################################
%% SOURCE RECONSTRUCTION %%

%% MEG preprocessing

cfg            = [];
cfg.dataset    = meg_file;
cfg.continuous = 'yes';
cfg.channel    = {'MEG'};
data           = ft_preprocessing(cfg);

cfg         = [];
cfg.length  = 2;
cfg.overlap = 0.5;
data        = ft_redefinetrial(cfg, data);

cfg        = [];
cfg.demean = 'yes';
data       = ft_preprocessing(cfg, data);

cfg            = [];
cfg.resamplefs = 100;
cfg.detrend    = 'yes';
meg_preprocessed = ft_resampledata(cfg, data);


%% Compute the leadfield
cfg = [];
cfg.grid = sourcemodel;
cfg.headmodel = headmodel;
cfg.channel = {'MEG'};
lf = ft_prepare_leadfield(cfg, meg_preprocessed);

%% MEG frequency analysis [frequency band]
% compute sensors level Fourier spectra (for cross-spectral density computation)
% from: https://www.fieldtriptoolbox.org/tutorial/beamformer/
%       https://github.com/fieldtrip/fieldtrip/blob/release/ft_freqanalysis.m

cfg           = [];
cfg.method    = 'mtmfft';
cfg.output    = 'fourier';
cfg.keeptrials = 'yes';
cfg.tapsmofrq = 1;
cfg.foilim    = [8 13];
freq_alpha           = ft_freqanalysis(cfg, meg_preprocessed);

%% Source reconstruction
cfg                   = [];
cfg.frequency         = [8 13];
cfg.method            = 'pcc';
cfg.grid              = lf;
cfg.headmodel         = headmodel;
cfg.sourcemodel       = sourcemodel;
cfg.keeptrials        = 'yes';
cfg.pcc.lambda        = '10%';
cfg.pcc.projectnoise  = 'yes';
cfg.pcc.fixedori      = 'yes';
source = ft_sourceanalysis(cfg, freq_alpha);

%% Compute connectivity
% from: https://www.fieldtriptoolbox.org/tutorial/networkanalysis/

cfg         = [];
cfg.method  = 'coh';
cfg.complex = 'absimag';
source_conn = ft_connectivityanalysis(cfg, source);

figure;
imagesc(source_conn.cohspctrm);

%% Read the atlas
atlas = ft_read_atlas('fieldtrip-20220707\template\atlas\aal\ROI_MNI_V4.nii');

%% Source interpolation

cfg = [];
cfg.parameter = 'tissue';
cfg.interpmethod = 'nearest';
sourcemodel_interpolated = ft_sourceinterpolate(cfg, atlas, sourcemodel);

%% Source parcelation

sourcemodel_interpolated.pos = source_conn.pos; % otherwise the parcellation won't work

cfg = [];
% cfg.parcellation = 'parcellation'; % ft_sourceparcellate(cfg, source_conn, sourcemodel_interpolated) -> Error: Unrecognized field name "parcellation".
cfg.parameter    = 'cohspctrm';
% parc_conn = ft_sourceparcellate(cfg, source_conn, atlas); % -> Error: the units of the source and parcellation structure are not consistent, please use FT_SOURCEINTERPOLATE
% parc_conn = ft_sourceparcellate(cfg, sourcemodel_interpolated, source_conn); % -> Error: This function requires parcellation data as input, see ft_datatype_parcellation.
parc_conn = ft_sourceparcellate(cfg, source_conn, sourcemodel_interpolated);

figure;
imagesc(parc_conn.cohspctrm);




Kind regards,
Kinga
This message and any attachment are intended solely for the addressee and may contain confidential information. If you have received this message in error, please contact the sender and delete the email and attachment. Any views or opinions expressed by the author of this email do not necessarily reflect the views of the University of Nottingham. Email communications with the University of Nottingham may be monitored where permitted by law.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20240503/6dadf411/attachment.htm>