Source reconstruction of event-related fields using minimum-norm estimation

Introduction

In this tutorial you can find information about how to do source-analysis with minimum-norm estimate on the event-related fields (MEG) of a single subject. We will working on the dataset described in the preprocessing tutorials (Trigger-based trial selection, Event related averaging and planar gradient), and we will use also the anatomical images that belong to the same subject. We will repeat code to select the trials and preprocess the data as described in the Event related averaging and planar gradient tutorial. We assume that preprocessing and event-related averaging is already clear for the reader. To preprocess the anatomical data, we will use two other software packages (FreeSurfer and MNE Suite).

This tutorial will not show how to do group-averaging and statistics on the source-level. It will also not describe how to do source-localization of oscillatory activation. You can check the Localizing oscillatory sources using beamformer techniques tutorial if you are interested in this latest.

Background

In the Event related averaging and planar gradient tutorial time-locked averages of event related fields of three conditions have been computed and the Cluster-based permutation tests on event related fields tutorial showed that there was a significant difference among two conditions. The topographical distribution of the ERFs belonging to each conditions and ERFs belonging to those differences have been plotted. The aim of this tutorial is to calculate a distributed representation of the underlying neuronal activity that resulted in the brain activity observed at the sensor level.

To calculate distributed neuronal activation we will use the minimum-norm estimate. This approach is favored for analyzing evoked responses and for tracking the wide-spread activation over time. It is a distributed inverse solution that discretizes the source space into locations on the cortical surface or in the brain volume using a large number of equivalent current dipoles. It estimates the amplitude of all modeled source locations simultaneously and recovers a source distribution with minimum overall energy that produces data consistent with the measurement 1) 2). The reference for the implemented method is Dale et al. (2000).

Procedure

Figure 1. shows the bigger steps in the calculation of the minimum-norm estimate. It shows that the computation of the inverse solution is based on the outputs of two independent processing steps: the processing of the anatomical images that leads to the forward solution and the processing of the MEG data. Creating the source model requires to use two additional software packages, FreeSurfer and MNE Suite.

Figure 1. An overview of the bigger steps in the calculation of the minimum-norm estimate

Figure 1. An overview of the bigger steps in the calculation of the minimum-norm estimate

To compute the distributed neuronal activation using minimum-norm estimate we will perform the following steps:

Processing of anatomical data


The following will use the anatomical MRI belonging to Subject01. The file can be obtained from ftp://ftp.fieldtriptoolbox.org/pub/fieldtrip/tutorial/Subject01.zip. The functions described in this part of the tutorial are using toolboxes that are under the fieldtrip/external folder. You do not have to add these toolboxes yourself, but it is important that you set up your MATLAB path properly. You can read about how to set up your MATLAB path here.

addpath <path to fieldtrip directory>
ft_defaults

Preprocessing of the anatomical MRI, reslicing and coregistration information


The following figure shows the steps of the preprocessing and their relation to the rest of the anatomical processing (volume conduction model and source model).
Figure 2. Processing of anatomical data
Figure 2.


The anatomical preprocessing is done in MATLAB with FieldTrip. The goal of this step is to create a file with a T1w anatomical image that can be used for the creation of two 'geometric objects': a volume conduction model, and a cortical sheet based source model. Moreover, this image will be used to create coregistration information to the different coordinate systems involved. One annoying thing to be aware of, and to think about in advance, is the fact that typically different parts of the pipeline assume (or require) different conventions of coordinate systems. Specifically, geometric information (sensor locations) in MEG/EEG data is typically expressed in a coordinate system that is defined based on external anatomical landmarks, whereas software used for processing of structural data usually requires coordinates to be expressed according to brain-anatomy related landmarks, such as the anterior and posterior commissures. More information about coordinate systems can be found at the frequently asked question about coordinate systems. To make a long story short, for now it suffices to know that you're safe if you know how to convert back and forth between the different relevant coordinate systems. The least error-prone and most convenient way to do this, is to create a well-defined reference anatomical image, which will be created with the following steps:

  • read in the anatomical images into MATLAB with ft_read_mri
  • ensure that the coordinate system is according to MNI's RAS convention, which can be checked with ft_determine_coordsys and imposed with ft_volumerealign.
  • reslice the volume with ft_volumereslice in order to have a uniform thickness for each slice, and to have the axes of the coordinate system lined up with the 'voxel axes'.
  • save the coregistration matrix that defines the transformation from voxels to MNI-RAS.
  • save the resliced anatomy in FreeSurfer compatible format, using ft_volumewrite. This anatomical image will serve as starting point for the creation of the cortical sheet based source model.

Then we may need to also create coregistration information to the coordinate system used in the MEG/EEG data.

  • realign the resliced anatomical data to the MEG/EEG based coordinate system with ft_volumerealign.
  • save the coregistration matrix that defines the transormation from voxels to the MEG/EEG coordinate system.

Thus, the input of the preprocessing is the anatomical MRI. The output is a realigned and resliced anatomical image, as well as a set of transformation matrices.


1. Preprocessing of the anatomical MRI: read in MRI data


mri_unknown = ft_read_mri('Subject01.mri');

Since by default we don't know in which coordinate system the MRI is specified, we call it mri_unknown.

2. Preprocessing of the anatomical MRI: impose coordinate system according to MNI convention


In this example, we are using an MRI which has been processed to contain a transformation matrix that corresponds to the CTF convention. As outlined above, life will be much easier if the coordinate system is in accordance with the MNI convention. Therefore, it needs to be 'realigned'. In general, it may be not clear what coordinate system is attached to the anatomical image.

To find out about the coordinate system of your mri, you can use the following function to check it:

mri_unknown = ft_determine_coordsys(mri_unknown, 'interactive', 'yes');

If it worked well, you will see the coordinate system specified in the mri structure in the mri_unknown.coordsys field. If 'coordsys' is not 'spm', you will need to align your mri to the anatomical landmarks (anterior commissure, posterior commissure, and a point that defines the postive z-direction) with the ft_volumerealign function. The mnemonic 'spm' for the coordsys is used to indicate this MNI-based convention. Ft_volumerealign does not change the anatomical data, instead it creates a transformation matrix that aligns the anatomical data to the intended coordinate system.

cfg = [];
cfg.method = 'interactive';
cfg.coordsys = 'spm';
mri_spm    = ft_volumerealign(cfg, mri_unknown);

3. Preprocessing of the anatomical MRI: reslicing


This step reslices the anatomical volume in a way that each slice will be equally thick. We use 1 mm thick slices and we specify the dimension as 256X256X256, because this is the format which FreeSurfer works with.

cfg            = [];
cfg.resolution = 1;
cfg.dim        = [256 256 256];
mrirs          = ft_volumereslice(cfg, mri_spm);
transform_vox2spm = mrirs.transform;

For convenience, you can now save the transformation matrix.

save('Subject01_transform_vox2spm', 'transform_vox2spm');

4. Preprocessing of the anatomical MRI: save to disk


% save the resliced anatomy in a FreeSurfer compatible format
cfg             = [];
cfg.filename    = 'Subject01';
cfg.filetype    = 'mgz';
cfg.parameter   = 'anatomy';
ft_volumewrite(cfg, mrirs);

Importantly, the mgz-filetype can only be used on Linux and Mac platforms. When you are processing the anatomical information on one of these platforms it is OK to save as mgz (and useful too, because it compresses the files and uses less diskspace as a consequence). These files cannot be saved nor read on a Windows PC. If you use MATLAB on Windows, you can save the volume as a nifti file using cfg.filetype = 'nifti'. Subsequently, if needed, you can convert it to mgz using mri_convert with FreeSurfer.


5. Preprocessing of the anatomical MRI: impose coordinate system according to M/EEG convention


Since this example is concerned with CTF-MEG data, we will need to get coregistration information that expresses anatomical information in coordinates according to the CTF-convention. To this end, we do a second realignment, now using the resliced anatomical image.

cfg = [];
cfg.method = 'interactive';
cfg.coordsys = 'ctf';
mrirs_ctf    = ft_volumerealign(cfg, mrirs);
transform_vox2ctf = mrirs_ctf.transform;

For convenience, you can now save the transformation matrix.

save('Subject01_transform_vox2ctf', 'transform_vox2ctf');

The MATLAB-based preprocessing of the anatomical data is now finished. We created an .mgz files that can be used for the creation of a cortical-sheet based source model and a volume conduction model. Moreover, 2 coregistration matrices have been constructed that allow to switch between coordinate systems.

Source model, creation of cortical sheet with Freesurfer + downsampling


Source model: Introduction


We will use FreeSurfer to create a source model that is based on a description of the cortical sheet. Essentially, we will construct a triangulated cortical mesh, ideally consisting of a number of approximately equally sized triangles that form a topological sphere for each of the cerebral hemispheres. The latter property is required to create an inflated cortex and to do intersubject realignment. FreeSurfer generates meshes with > 100000 vertices per hemisphere, which is too much for a workable M/EEG source reconstruction. Therefore, we use the MNE-suite to downsample the triangulated meshes. This step serves the purpose of retaining a topologically correct description of the surface, and keeping the variance in triangle size low. In contrast, MATLAB's reducepatch function breaks the topology and leads to a bigger variance in triangle size.
The creation process of the source-space can be divided into 4 stages (after the preprocessing steps):

  1. Volumetric processing in FreeSurfer.
  2. Surface based processing in FreeSurfer.
  3. Creation of the mesh using MNE-suite.
  4. Coregistration of the source space to the sensor-based coordinate system with FieldTrip.

The volumetric and surface based processing (first and the second steps) can be together 10 hours long. These steps will run on the computer from themselves. There is only one checkpoint between the volume and the surface based processing when an intermediate result can be checked interactively.

The input of the creation process of the source-space are mgz files that were created in the preprocessing. The output is the source model that is a MATLAB structure called 'sourcespace' in this tutorial.

The instructions about how to install and run FreeSurfer and MNE Suite are aimed at users at the Center of Neuroimaging of the Donders Institute and at the MPI for Psycholinguistics in Nijmegen.

1. Source model: Volumetric processing in FreeSurfer


FreeSurfer's anatomical processing pipeline consists of a series of automated steps, which essentially consist of:

  1. processing steps on a volumetric anatomical MRI (image intensity normalization, co-registration with Talairach space, skull stripping, automatic segmentation of sub-cortical structures, and finally segmentation)
  2. extraction of the cortical mesh
  3. processing of the surface meshes (smoothing, topology fixing, inflation, co-registration with a spherical template)

Here is a link to the different processing steps. Although the FreeSurfer procedure can be invoked using only a few FreeSurfer commands, below we will describe the (sub)commands that will achieve the same. These commands sequentially generate a series of files (volumetric, surface and transformation matrices). Each of the output files serves as input to the sequential analysis steps. A table of file dependencies can be found here.
There are a few analysis steps in FreeSurfer which are not guaranteed to give a nice result, and require some user interaction to get it right. Moreover, FreeSurfer can be quite picky with respect to the exact format of the MRI-volumes. One step which in our experience is notorious for not being very robust is automatic skull-stripping. Therefore, we advocate a hybrid approach that uses SPM for an initial segmentation of the anatomical MRI during the preprocessing. With this segmentation, we can create a skull-stripped image, which is a prerequisite for a correct segmentation in FreeSurfer. Although this approach may seem a bit convoluted (you may rightfully ask why we need to redo the segmentation in FreeSurfer if we already did it in SPM), the interdependencies between different files generated along the FreeSurfer pipeline make tapping into this pipeline at a random point quite complicated. For this reason a large part of the volumetric processing in FreeSurfer needs to be done as well.

To create a skullstripped anatomical image, you can do the following. We assume that you have executed all steps that are described in the section about preprocessing of the anatomical MRI, and that you have a file that contains the resliced anatomical image, expressed in the MNI-RAS coordinate system.

mri = ft_read_mri('Subject01.mgz');
mri.coordsys = 'spm';

cfg = [];
cfg.output = 'brain';
seg = ft_volumesegment(cfg, mri);
mri.anatomy = mri.anatomy.*double(seg.brain);

cfg             = [];
cfg.filename    = 'Subject01masked';
cfg.filetype    = 'mgz';
cfg.parameter   = 'anatomy';
ft_volumewrite(cfg, mri);


In order to be able to use FreeSurfer, you need to have a working installation of the package. It can be downloaded from here. If you are working at the Center of Neuroimaging of the Donders Institute you can find more versions of FreeSurfer under the /opt/FreeSurferXXX directories. (If you are working at the MPI for Psycholinguistics, you should install the software yourself in your directory.) We recommend to use FreeSurfer 5.3. You can run the commands just copying and pasting them into the terminal window of the Linux system (from where you used also MATLAB).

To get started, you need to set up your environment variables. Please pay close attention to the spaces in the following commands, or the lack thereof.

export FREESURFER_HOME=<path to FreeSurfer>
export SUBJECTS_DIR=<Subject directory>

SUBJECTS_DIR is the directory where you will store all the FreeSurfer-processed anatomical data of all your subjects. Then, type this command to set up FreeSurfer:

source $FREESURFER_HOME/SetUpFreeSurfer.sh

The following populates an empty subject-specific directory with subdirectories:

mksubjdirs $SUBJECTS_DIR/Subject01

Now, we are ready to start using FreeSurfer. As a first step in the volumetric pipeline, we have to 'convert' the anatomical MRI once more, but now using a FreeSurfer command. You start by making a new folder in the subject directory called “mri” into which you will copy both the masked and the original mgz files you created in the previous preprocessing steps in FieldTrip. All subsequent FreeSurfer commands will be called from the “mri” directory.

cp Subject01masked.mgz $SUBJECTS_DIR/Subject01/mri/Subject01masked.mgz
cp Subject01.mgz       $SUBJECTS_DIR/Subject01/mri/Subject01.mgz
cd $SUBJECTS_DIR/Subject01/mri

mri_convert -c -oc 0 0 0 Subject01masked.mgz brainmask.mgz
mri_convert -c -oc 0 0 0 Subject01.mgz       orig.mgz

recon-all -talairach      -subjid Subject01
recon-all -nuintensitycor -subjid Subject01
recon-all -normalization  -subjid Subject01
recon-all -gcareg         -subjid Subject01
recon-all -canorm         -subjid Subject01
recon-all -careg          -subjid Subject01
recon-all -careginv       -subjid Subject01
recon-all -calabel        -subjid Subject01
recon-all -normalization2 -subjid Subject01
recon-all -maskbfs        -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:

cd <Subject directory>/Subject01/mri
mri = ft_read_mri('filled.mgz');

cfg = [];
cfg.interactive = 'yes';
ft_sourceplot(cfg, mri);

Figure 3. Filled mgz


Figure 3. Filled mgz created by FreeSurfer. The two hemispheres have different colors (white and grey), cerebellum is not included.

2. Source model: Surface based processing in FreeSurfer


The surface construction is done by the following sequence of commands (from the Subject01/mri directory):

recon-all -fill       -subjid Subject01
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
recon-all -white      -subjid Subject01
recon-all -finalsurfs -subjid Subject01
recon-all -smooth2    -subjid Subject01
recon-all -inflate2   -subjid Subject01

# then use a shortcut command to do the rest, but we need the rawavg.mgz file to exist
cp $SUBJECTS_DIR/Subject01/mri/Subject01.mgz $SUBJECTS_DIR/Subject01/mri/rawavg.mgz
recon-all -autorecon3 -subjid Subject01

After these steps (which may take quite a while) you end up with a bunch of files in the Subject01/surf/ directory. We are going to use lh.white and rh.white to create the source space in the next step.

3. Source model: Creation of the mesh using MNE Suite


Just like with FreeSurfer, we have to first take care that MNE-suite is installed, and that some environmental variables are correctly specified. If you are working in Nijmegen at the DCCN, you can find the MNE suite under the /opt/mne directory. At the MPI, the MNE suite is installed under the /mnt/data1/mne directory.

export MNE_ROOT=<path to MNE>
source $MNE_ROOT/bin/mne_setup.sh

export SUBJECTS_DIR=<Subject directory>
export SUBJECT=Subject01

Now we can create the source space

mne_setup_source_space --ico -6

This step creates a bunch of files in <Subject directory>/Subject01/bem/, containing different representations of the source space. In subsequent steps, FieldTrip will use the Subject01-oct-6-src.fif file. We can already have a look in MATLAB at how the source space looks.

sourcespace = ft_read_headshape('Subject01-oct-6-src.fif', 'format', 'mne_source');

figure
ft_plot_mesh(sourcespace);

Figure 4. Source-space downsampled
Figure 4. The source-space downsampled by MNE Suite

4. Source model: Co-registration of the source space to the sensor-based head coordinate system


We have the source locations co-registered to the MNI coordinate system, so now we need to co-register the source space to the sensor-array (i.e., we have to express the positions of the sources in the same coordinate system as the sensors). For this, we will use the transformation matrices computed in earlier in this tutorial. Specifically, using the resliced anatomical data in the mrirs-structure, we obtained a set of 2 transformation matrices, that describe the mapping of anatomical volumetric voxel indices into the sensor-based coordinate system transform_vox2ctf and into the MNI coordinate system transform_vox2spm. These two matrices can be combined in the following way, to yield a transformation matrix that transforms from MNI coordinates to sensor-based coordinates:

T = transform_vox2ctf/transform_vox2spm;

Now we can use this transformation matrix, to get it in the correct coordinate system.

% go to the Subject01/bem directory
sourcespace = ft_read_headshape('Subject01-oct-6-src.fif', 'format', 'mne_source');
sourcespace = ft_convert_units(sourcespace, 'mm');
sourcespace = ft_transform_geometry(T, sourcespace);

save sourcespace sourcespace

Volume conduction model


The volume conduction model describes the geometry and the electrical (conductive) properties of the head. The volume conduction model (or headmodel) requires a geometrical description of the head. In this tutorial we have an anatomical MRI of the subject from which we can construct the head model. However, it is important that the anatomical MRI and the sensor positions are expressed in the same coordinate system. The MEG sensor positions are always defined relative to the fiducial coils. If we want to create a volume conduction model, the anatomical data must also be expressed relative to these points (i.e., in the CTF coordinate system). Here, you can read more about the coordinate systems. And here, you can read more about the headmodel.

We create the volume conduction model from Subject01.mgz file. If we make sure to use the sensor-based coordinate system in this anatomical image, we should be alright.

mri = ft_read_mri('Subject01.mgz');
load('Subject01_transform_vox2ctf.mat');

mri.transform = transform_vox2ctf; % this is only allowed, if the transformation matrix was constructed in the exact same anatomical image as the one in the 'Subject01.mgz' file.
mri.coordsys  = 'ctf';

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

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

It is useful to check if the resulting sourcespace and the volume conductor are aligned. To plot, you can use this code:

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

If they are not aligned, it may be because vol is not expressed in CTF coordinates. You can check using ft_determine_coordsys.

Figure 5. Source-space with volume conductor
Figure 5. The final version of the source-space aligned and plotted together with the volume conductor

Processing of functional data


The following will use the MEG data belonging to Subject01. The file can be obtained from ftp://ftp.fieldtriptoolbox.org/pub/fieldtrip/tutorial/Subject01.zip. For both preprocessing and averaging, we will follow the steps that have been written in the Event related averaging and planar gradient tutorial. We will use trials belonging to two conditions (FC and FIC) and we will calculate their difference.

Preprocessing of MEG data

Reading the FC data

The ft_definetrial and ft_preprocessing functions require the original MEG dataset, which is available from ftp://ftp.fieldtriptoolbox.org/pub/fieldtrip/tutorial/Subject01.zip.

% find the interesting segments of data
cfg = [];                                           % empty configuration
cfg.dataset                 = 'Subject01.ds';       % name of CTF dataset  
cfg.trialdef.eventtype      = 'backpanel trigger';
cfg.trialdef.prestim        = 1;
cfg.trialdef.poststim       = 2;
cfg.trialdef.eventvalue     = 9;                    % trigger value for fully congruent (FC)
cfg = ft_definetrial(cfg);            

% remove the trials that have artifacts from the trl
cfg.trl([2, 3, 4, 30, 39, 40, 41, 45, 46, 47, 51, 53, 59, 77, 85],:) = []; 

% preprocess the data
cfg.channel    = {'MEG', '-MLP31', '-MLO12'};       % read all MEG channels except MLP31 and MLO12
cfg.demean     = 'yes';
cfg.baselinewindow  = [-0.2 0];
cfg.lpfilter   = 'yes';                              % apply lowpass filter
cfg.lpfreq     = 35;                                 % lowpass at 35 Hz.

dataFC_LP = ft_preprocessing(cfg);                      

These data have been cleaned from artifacts by removing several trials and two sensors; see the visual artifact rejection tutorial.

Subsequently you can save the data to disk.

save dataFC_LP dataFC_LP

Reading the FIC data

The ft_definetrial and ft_preprocessing functions require the original MEG dataset, which is available from ftp://ftp.fieldtriptoolbox.org/pub/fieldtrip/tutorial/Subject01.zip.

% find the interesting segments of data
cfg = [];                                           % empty configuration
cfg.dataset                 = 'Subject01.ds';       % name of CTF dataset  
cfg.trialdef.eventtype      = 'backpanel trigger';
cfg.trialdef.prestim        = 1;
cfg.trialdef.poststim       = 2;
cfg.trialdef.eventvalue     = 3;                    % trigger value for fully incongruent (FIC)
cfg = ft_definetrial(cfg);            

% remove the trials that have artifacts from the trl
cfg.trl([15, 36, 39, 42, 43, 49, 50, 81, 82, 84],:) = []; 

% preprocess the data
cfg.channel    = {'MEG', '-MLP31', '-MLO12'};        % read all MEG channels except MLP31 and MLO12
cfg.demean     = 'yes';
cfg.baselinewindow  = [-0.2 0];
cfg.lpfilter   = 'yes';                              % apply lowpass filter
cfg.lpfreq     = 35;                                 % lowpass at 35 Hz.

dataFIC_LP = ft_preprocessing(cfg);                      

These data have been cleaned from artifacts by removing several trials and two sensors; see the visual artifact rejection tutorial.

Subsequently you can save the data to disk.

save dataFIC_LP dataFIC_LP

Averaging and noise-covariance estimation

The function ft_timelockanalysis makes averages of all the trials in a data structure and also estimates the noise-covariance. For a correct noise-covariance estimation it is important that you used the cfg.demean = 'yes' option when the function ft_preprocessing was applied.

The trials belonging to one condition will now be averaged with the onset of the stimulus time aligned to the zero-time point (the onset of the last word in the sentence). This is done with the function ft_timelockanalysis. The input to this procedure is the dataFC_LP structure generated by ft_preprocessing. At the same time, we need to compute the noise-covariance matrix, therefore cfg.covariance = 'yes' has to be specified as well as the time window where the noise-covariance will be estimated. Here, we use the baseline where there is no signal of interest yet.

  load dataFC_LP;
  load dataFIC_LP;
  cfg = [];
  cfg.covariance = 'yes';
  cfg.covariancewindow = [-inf 0]; %it will calculate the covariance matrix 
                                   % on the timepoints that are  
                                   % before the zero-time point in the trials
  tlckFC = ft_timelockanalysis(cfg, dataFC_LP);
  tlckFIC = ft_timelockanalysis(cfg, dataFIC_LP);
  save tlck tlckFC tlckFIC;

Forward solution

The source space, the volume conduction model and the position of the sensors are necessary inputs for creating the leadfield (forward solution) with the ft_prepare_leadfield function. The sensor positions are contained in the grad field of the averaged data. However, the grad field contains the positions of all channels, therefore, the used channels have to be also specified.

load tlck;
load sourcespace;
load vol;

cfg = [];
cfg.grad = tlckFC.grad;                      % sensor positions
cfg.channel = {'MEG', '-MLP31', '-MLO12'};   % 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;                               % volume conduction model
leadfield = ft_prepare_leadfield(cfg);

save leadfield leadfield;

Inverse solution


The ft_sourceanalysis function calculates the inverse solution. The method used (minimum-norm estimation) has to be specified with the cfg.method option. The averaged functional data, the forward solution (the output of the ft_prepare_leadfield function), the volume conduction model (in this case, the output of the ft_prepare_headmodel function) and the noise-covariance matrix (the cov field of the output of the ft_timelockanalysis function) have to be provided.

The lambda value is a scaling factor that is responsible for scaling the noise-covariance matrix. If it is zero the noise-covariance estimation will be not taken into account during the computation of the inverse solution. Noise-covariance is estimated in each trial separately and then averaged, while the functional data (of which we calculate the source-analysis) is simply averaged across all the trials. Therefore, the higher the number of trials the lower the noise is in the averaged, functional data, but the number trials is not reducing the noise in the noise-covariance estimation. This is the reason while it is useful to use a scaling factor for the noise-covariance matrix if we want to estimate more realistically the amount of noise.

You do not have to specify of the noise-covariance matrix separatly, because it is in the tlckFC.cov and in the tlckFIC.cov fields, and ft_sourceanalysis will take it into account automatically.

load tlck;
load leadfield;
load vol;

cfg        = [];
cfg.method = 'mne';
cfg.grid   = leadfield;
cfg.vol    = vol;
cfg.mne.prewhiten = 'yes';
cfg.mne.lambda    = 3;
cfg.mne.scalesourcecov = 'yes';
sourceFC  = ft_sourceanalysis(cfg,tlckFC);
sourceFIC = ft_sourceanalysis(cfg, tlckFIC);

save source sourceFC sourceFIC;

Visualization

You can plot the inverse solution onto the source-space at a specific time-point with the ft_plot_mesh function.

load source;
load sourcespace;

bnd.pnt = sourcespace.pnt;
bnd.tri = sourcespace.tri;
m=sourceFIC.avg.pow(:,450); % plotting the result at the 450th time-point that is 
                         % 500 ms after the zero time-point
ft_plot_mesh(bnd, 'vertexcolor', m);

Figure 6. The source reconstruction at 500 ms
Figure 6. The result of the source-reconstruction of the FIC condition plotted onto the source-space at 500 ms after the 0 time-point

But we would like to know where the difference between the conditions can be localized. Therefore, we calculate the difference of the two conditions, and we use ft_sourcemovie to visualize the results.

cfg = [];
cfg.projectmom = 'yes';
sdFC = ft_sourcedescriptives(cfg,sourceFC);
sdFIC = ft_sourcedescriptives(cfg, sourceFIC);

sdDIFF = sdFIC;
sdDIFF.avg.pow = sdFIC.avg.pow - sdFC.avg.pow;
sdDIFF.tri = sourcespace.tri;

save sd sdFC sdFIC sdDIFF;

cfg = [];
cfg.mask = 'avg.pow';
ft_sourcemovie(cfg,sdDIFF);

Figure 7. One frame from ft_sourcemovie
Figure 7. One frame from the movie that shows the differences of the two source-reconstructions

Summary and further readings

In this tutorial we showed how to do MNE source reconstruction method on a single subject data. We compared the averaged ERF in two conditions and we reconstructed the sources and we calculated the difference of the two source-reconstruction. We showed also how you can visualize the results.

Functions and tutorial pages that show how to average, and how to analyze statistically source-reconstructions across subjects or how to compare those to a template brain are still under development.

FAQs:

Example scripts:

1) Ou, W., Hamalainen, M., Golland, P., 2008, A Distributed Spatio-temporal EEG/MEG Inverse Solver
2) Jensen, O., Hesse, C., 2010, Estimating distributed representation of evoked responses and oscillatory brain activity, In: MEG: An Introduction to Methods, ed. by Hansen, P., Kringelbach, M., Salmelin, R., doi:10.1093/acprof:oso/9780195307238.001.0001
tutorial/minimumnormestimate.txt · Last modified: 2016/10/24 13:46 by 131.174.45.112

You are here: starttutorialminimumnormestimate
CC Attribution-Share Alike 3.0 Unported
www.chimeric.de Valid CSS Driven by DokuWiki do yourself a favour and use a real browser - get firefox!! Recent changes RSS feed Valid XHTML 1.0
This DokuWiki features an Anymorphic Webdesign theme, customised by Eelke Spaak and Stephen Whitmarsh.
Mobile Analytics Website Security Test