Tutorial 22: Source estimation [TODO]

Authors: Francois Tadel, Elizabeth Bock, Rey R Ramirez, John C Mosher, Richard Leahy, Sylvain Baillet

You have in your database a forward model that explains how the cortical sources determine the values on the sensors. This is useful for simulations, but what we need is to build the inverse information: how to estimate the sources when we have the recordings. This tutorials introduces the tools available in Brainstorm for solving this inverse problem.

WARNING: The new interface presented here do not include all the options yet. The sLORETA normalization or the mixed head models are not supported: to use them, use the old interface (menu "Compute sources" instead of "Compute sources 2016", described in the ?old tutorials).

Ill-posed problem

Our goal is to estimate the activity of the thousands of dipoles described by our forward model. However we only have a few hundred variables in input (the number of sensors). This inverse problem is ill-posed, meaning there is an infinite number of combinations of source activity patterns that can generate exactly the same sensor topography. Inverting the forward model directly is impossible, unless we add some strong priors to our model.

Wikipedia says: "Inverse problems are some of the most important and well-studied mathematical problems in science and mathematics because they tell us about parameters that we cannot directly observe. They have wide application in optics, radar, acoustics, communication theory, signal processing, medical imaging, computer vision, geophysics, oceanography, astronomy, remote sensing, natural language processing, machine learning, nondestructive testing, and many other fields."

Many solutions have been proposed in the literature, based on different assumptions on the way the brain works and depending on the amount of information we already have on the effects we are studying. Among the hundreds of methods available, two classes of inverse models have been widely used in MEG/EEG source imaging in the past years: minimum-norm solutions and beamformers.

Both approaches have the advantage of being linear: the activity of the sources is a linear recombination of the MEG/EEG recordings. It is possible to solve the inverse problem independently of the recordings, making the data manipulation a lot easier and faster.

Both are available in Brainstorm, so you can use the one that is most adapted to your recordings or to your own personal expertise. Only the minimum norm estimates will be described in this tutorial, but the other solutions work exactly the same way.

Source estimation options [TODO]

Before we start estimating the sources for the recordings available in our database, let's start with an overview of the options available. The screen capture below represents the basic options for the minimum norm estimates. The options for the other methods will be described in advanced tutorials.


Method [TODO]

Measure [TODO]

The minimum norm estimates produce a measure of the current density flowing at the surface of the cortex. To visualize these results and compare them between subjects, we can normalize the MNE values to get a standardized level of activation with respect to the noise or baseline level (dSPM, sLORETA, etc).

Source orientation

Computing sources for an average

Display: Cortex surface

Why does it look so noisy?

The source maps look very noisy and discontinuous, they show a lot of disconnected patches. This is due to the orientation constraint we imposed on the dipoles orientations. Each value on the cortex has to be interpreted as a vector, oriented perpendicular to the surface. Because of the brain circumvolutions, all the sources have different orientations, two adjacent sources have very little chance to have the same orientation in this model, therefore the minimum norm method may attribute completely different values to them. This causes all these gaps we see here.

Visually, you should not always interpret disconnected colored patches as independent sources. You cannot expect a very spatial resolution with this technique (~5-10mm). Most of the time, a cluster of disconnected source patches in the same neighborhood that show the same evolution in time can be interpreted as "there is some significant activity around here, but we don't know where exactly".

To get more continuous maps for visualization or publication purposes, you can either smooth the values explicitly on the surface (process "Sources > Spatial smoothing") or use unconstrained source models.

For data exploration, this is a good enough representation of the brain activity, mostly because it is fast and efficient. You can get a better feeling of the underlying brain activity patterns by making short interactive movies: click on the figure, then hold the left or right arrows of your keyboard.

Activity patterns will also look sharper when we compute normalized measures (later in this tutorial). In most of the screen captures in this following sections, the contrast of the figures has been enhanced for illustration purposes. Don't worry if it looks a lot less colorful on your screen.

Display: MRI Viewer

Display: MRI 3D

Sign of constrained maps

You should pay attention to the sign of the current amplitudes that are given by the minimum norm method: they can be positive or negative and they oscillate around zero. Display the sources on the surface, set the amplitude threshold to 0%, then configure the colormap to show relative values (uncheck the "Absolute values" option), you would see those typical stripes of positive and negative values around the sulci. Double-click on the colorbar after testing this to reset the colormap.

This pattern is due to the orientation constraint imposed on the dipoles. On both sides of a sulcus, we have defined dipoles that are very close to each other, but with opposite orientations. If we have a pattern of activity on one side of a suclus that can be assimilated to an electric dipole (green arrow), the minimum norm model will try to explain it with the dipoles that are available in the head model (red and blue arrows). Because of the dipoles orientations, it translates into positive values (red arrows) on one side of the sulcus and negative on the other side (blue arrows).

When displaying the cortical maps at one time point, we are usually not interested by the sign of the minimum norm values but rather by their amplitude. This is why we always display them by default with the colormap option "absolute values" selected.

However, we cannot simply discard the sign of these values because we need them for other types of analysis, typically time-frequency decompositions and connectivity analysis. For estimating frequency measures on the source maps, we need to keep the oscillations around zero.

Unconstrained orientations

In the cases where the orientation constraint imposed on the dipoles orientations looks too strong, it is possible to relax it partially (option "loose constraints") or completely (option "unconstrained"). This is typically something to consider when using a MNI template instead of the subject's anatomy, or when studying deeper or non-cortical brain regions for which the normal to the FreeSurfer cortex surface is unlikely to match any physiological reality.

In terms of data representation, the option "unconstrained" and "loose constraints" are very similar. Instead of using one dipole at each cortical location, a base of three orthogonal dipoles is used.
Here we will only illustrate the fully unconstrained case.

Source map normalization

The current density values returned by the minimum norm method have a few problems:

Normalizing the current density maps with respect to a reference level (estimated from noise recordings, pre-stimulus baseline or resting state recordings) can help with all these issues at the same time. Some normalizations can be computed independently from the recordings, and added to the linear inverse operator (dSPM or sLORETA). Another way of proceeding is to divide the current density maps by the standard deviation estimated over a baseline (Z-score).

The normalization options do not change the temporal dynamics of your results, they are just different ways for looking at the same minimum norm maps. If you look at the time series associated with one given source, it would be exactly the same for all the normalizations, except for a scaling factor. Only the relative weights change between the sources, and these weights do not change over time.



Typical recommendations

Delete your experiments

Computing sources for single trials

Because the minimum norm model is linear, we can compute an inverse model independently from the recordings and apply it on the recordings when needed. We will now illustrate how to compute a shared inverse model for all the imported epochs.

Averaging in source space

Computing the average

Visualization filters

Z-score normalization


Averaging normalized values

Averaging normalized source maps within a single subject requires more attention than averaging current density maps. The amplitude of the normalized measures increase with the SNR of the signal, the higher the SNR the higher the normalized score. For instance, the average of the dSPM for the single trials is lower than the dSPM of the averaged trials.





Display: Contact sheets and movies

A good way to represent what is happening in time is to generate contact sheets or videos. Right-click on any figure and go to the menu Snapshot to check out all the possible options. For a nicer result, take some time to adjust the size of the figure, the amplitude threshold and the colormap options (hiding the colorbar can be a good option for contact sheets).

A time stamp is added to the captured figure. The size of the text font is fixed, so if you want it to be readable in the contact sheet, you should make you figure very small before starting the capture.


Model evaluation

One way to evaluate the accuracy of the source reconstruction if to use simulate recordings using the estimate source maps. This is done simply by multiplying the source time series with the forward model:
MEG_simulated [Nmeg x Ntime] = Forward_model [Nmeg x Nsources] * MN_sources [Nsources x Ntime]
Then you can compare the original MEG recordings with the simulated ones visually. More formally, you can compute an error measure from the residuals (recordings - simulated).

To simulate MEG recordings from a minimum norm source model, right-click on the source file, then select the menu "Model evaluation > Simulate recordings".

Open side-by-side the original and simulated MEG recordings for the same condition:


Advanced options [TODO]

Right-click on the deviant average for Run#01 > Compute sources [2016].
Click on the button [Show details] to bring up all the advanced minimum norm options.

Depth weighting [TODO]

The minimum-norm estimates have a bias towards superficial currents. This tendency can be alleviated by adjusting these parameters. To understand how to set these parameters, please refer to the MNE manual (options --depth, --weightexp and --weightlimit).

Noise covariance regularization [TODO]

Explain the influence of this parameter.

Regularization parameter [TODO]

Explain the influence of this parameter.

Output mode [TODO]

LCMV beamformer [TODO]

Explain the method briefly.

Equations [TODO]



On the hard drive [TODO]

Constrained shared kernel

Right-click on a shared inverse file in the database explorer > File > View file contents.

Structure of the source files: results_*.mat

Mandatory fields:

Optional fields:

Full source maps

In Intra-subject, right-click on one of the normalized averages > File > View file contents.

This file has the same structure as a shared inverse kernel, with the following differences:

Filename tags

Useful functions

References [TODO]

Additional documentation

Feedback: Comments, bug reports, suggestions, questions
Email address (if you expect an answer):

Tutorials/SourceEstimation (last edited 2016-04-11 17:45:38 by FrancoisTadel)