SPM Source reconstruction demo
During this demo we will estimate the current sources on the cortical surface that gave rise to the electromagnetic activity we measured with MEG and EEG outside the head. To do so, we first need to build a head model, based on the individual MRI, coregister the head model and our M/EEG data and then build a forward model which describes how sources in the brain appear on the different M/EEG sensors. We will then apply two different inversion approaches, both in a Bayesian framework: the Multiple Sparse Priors and the Minimum Norm solution. We will compare the model evidences for these distributed imaging approaches and source localise the evoked oscillatory activity at a particular time-frequency window.
We start by opening SPM1
-
First go to the source reconstruction demo folder: cd C:/workshop/spm-source-demo.
-
Type âspm eegâ in the command line to start SPM.
Create Head Model and Forward Model
We will use the subjectâs structural MRI to define meshes describing the cortex, skull and scalp of our subject. Due to the different conductivities (for example between the skull and the scalp) the geometry of the head is important to derive realistic physical forward models that predict how cortical current sources map onto the M/EEG sensors.
- In the SPM window click on the batch button. In the batch window select âSPM â M/EEG - Source reconstruction - Head model specificationâ
- Double-click on âM/EEG datasetsâ and select âPapMcbdspmeeg_run_01_sssâ from the Sensor-level statistics demo as input file. This file contains artifact corrected single trial data for EEG, MEG magnetometers, planar gradiometers and combined gradiometers.
- Next, set the inversion index to â3â. This index allows you to track different types of forward models and inverse solutions and compare their log model evidences. We choose â3â here as we will later use some pre-calculated inversion models filed under the inversion indices 1 and 2.
- Additional comments relating to each index can be inserted if âcommentsâ is selected. We leave this field empty for now.
- In the field âmeshesâ first select âmesh sourceâ. From here specify âIndividual structural imageâ and select the NIfTI file âmprage.niiâ. The mesh resolution can be kept as normal (approximately 4000 vertices per hemisphere). This step will generate scalp, skull and cortical meshes by warping meshes from a template brain with the inverse spatial normalisation of the subjectâs brain.
- To coregister the MRI and M/EEG data, please select âspecify coregistration parametersâ. We first need to specify the fiducials, with coordinates in the space of the MRI image selected. Here we will define âNasionâ and the left and right pre-auricular points âLPAâ and âRPAâ.
- Set the first M/EEG fiducial label to âNasionâ.
- For the Nasion, highlight âHow to specifyâ, select âtype MRI coordinatesâ and enter â[2 81 -4]â
- Then specify the second fiducial label as âLPAâ and enter â[-77 -14 -23 ]â for MRI coordinates; the final fiducial label is âRPAâ, with MRI coordinates â[83 -20 -25]â.
- Select ânoâ for the âuse headshape pointsâ option.
- Finally, for the forward model itself, we highlight âEEG head modelâ, and specify this as âEEG BEMâ; we select âMEG head modelâ and specify this as âSingle Shellâ.
You can now run this batch file by pressing the green ârunâ button. Note that at this stage the model parameters are saved, but the lead fields themselves are not estimated until inversion.
Review the head model and forward model
- You can review the head model and the forward model by pressing the â3D Source Reconstructionâ button within the SPM window.
- This will create a new window; select âLoadâ and choose the âPapMcbdspmeeg_run_01_sss.matâ file. You can find the forward model we just computed under the inversion index â3â (use the ânextâ and âpreviousâ buttons to navigate between different inversion indices).
- On the left hand side of the âsource localisationâ window, select the âdisplayâ button below the âMRIâ button. This will bring up the scalp (orange), inner and outer skull (red) and cortical (blue) meshes of the subjectâs brain. Note that the fiducials are shown by cyan discs.
- Next, select the âdisplayâ button beneath âCo-registerâ and then select âMEGâ when asked what to display. The graphics window should the sensor locations in green discs, the digitized headpoints in small red dots, the fiducials in the MEG data as purple diamonds, and the MRI fiducials as cyan discs again. The overlap between MRI and MEG fiducials indicates how well the data have been coregistered.
- Finally, select âDisplayâ, located beneath âForward Modelâ, and select again âMEGâ. You should see an image displaying the MEG sensor locations relative to the cortical mesh and the single shell fitted to the inner skull mesh.
Model Inversion
We will next compute the inverse solutions and compare two distributed source reconstruction approaches. Since the inverse problem is ill-posed, prior information must be included to give a unique solution. We will first apply the Multiple Sparse Priors approach; this corresponds to a sparse prior on the sources, namely that only a few sources are activ
- Go back to the Batch Editor, under âfileâ open a new batch and select âSPM - M/EEG - Source reconstruction â Source Inversionâ
- Select the same input file âPapMcbdspmeeg_run_01_sss.matâ, and set the inversion index to â1â. We will make use of some pre-calculated inversion models created earlier for indices 1 and 2.
- We will invert the data for Famous, Unfamiliar and Scrambled faces: highlight âwhat conditions to includeâ and select âAllâ.
- Next highlight inversion parameters, choose âcustomâ and keep the inversion type as âGSâ. Greedy Search (GS) is one of several fitting algorithms for optimising the MSP approach; we choose GS here because it is quickest and works well for these data.
- Then enter the time window of interest as â[-100 800]â
- We want to consider all frequencies; please set the frequency window of interest to â[0 256]â.
- Select âyesâ for the âPST Hanning windowâ. In this example we do not provide âSource priorsâ and donât want to âRestrict solutionsâ; these fields remain unchanged. You can find an example on how to include fMRI-based location priors into source reconstruction in the SPM12 manual chapter on multimodal and multisubject data fusion.
- As final step use âSelect Modalitiesâ and choose âMEGâ to invert only the 102 magnetometers. You can also select several sensor types and fuse EEG and MEG data. Note that we cannot invert the combined planar gradients, because we do not have a physical forward model for those.
The second type of inversion we will examine corresponds to a L2-minimum norm (MNM) solution, i.e, fitting the data at the same time as minimising the total energy of the sources. In SPM, this is called âIIDâ because it corresponds to assuming that the prior probability of each source being active is independent and identically distributed.
- Go back to the batch editor, add another âM/EEG - Source reconstruction â Source Inversionâ module, and select the same input files as before (âPapMcbdspmeeg_run_01_sss.matâ), but this time set the inversion index to 2. We use the same forward model as for the âMSPâ source reconstruction.
- Set the inversion parameters to âcustomâ again, but this time with the inversion type âIIDâ. The remaining parameters should be made to match the Multiple Sparse Priors inversion approach.
Time-frequency contrasts
With the previous steps, we invert data of whole trials from -100 to +800ms across all frequencies. Next we want to localise the evoked activity around the N/M170 face component by averaging power across a suitable time-frequency window
- Select âM/EEG - Source reconstruction â Inversion Resultsâ.
- Again, select âPapMcbdspmeeg_run_01_sss.matâ as input file; set the inversion index to 1. We will first investigate averaged source power across the time-frequency window of interest for the MSP solution.
- Set the time window of interest to â[100 250]â and the frequency window of interest to â[10 20]â.
- For the contrast type, please select âevokedâ from the current item window, and the output space as âMNIâ.
- We can now write the source power image as a surface-based GIFTI âMeshâ by high-lighting the option âMeshâ.
- Finally, we smooth across the cortical surface to render our data more suitable for RFT. Keep the default cortical smoothing of 8.
- Then replicate this module to produce a second âinversion resultsâ module by right-clicking on the module and selecting âreplicateâ.
- Here, change the inversion index from â1â to â2â to write out the time-frequency contrast for the MNM (IID) distributed solution as well.
Generating batch scripts
You can use the batch framework to construct a processing pipeline across a group of subjects or across different source reconstruction approaches. A brief example
- Save the batch by going to âFileâ in the Batch editor and select âsave Batch and Scriptâ. As file name write âbatch_localise_invâ. This will result in two files: a batch file âbatch_localise_inv_job.mâ and a Matlab script âbatch_localise_inv.mâ which runs the batch file. Take a look at the job-file by entering âopen batch_localise_inv_job.mâ in the MATLAB command window.
- In batch_localise_inv.m replace nrun = X; with nrun = 1 and save the file. We could now run this Matlab script with the same results as if we had pressed the green âRunâ button in the Batch Editor.
- To make things a bit more interesting, we go back to the batch and right-click on the M/EEG datasets field of the first source inversion module and select âClear Valueâ. Do the same for âTime window of interestâ and save batch and script under âbatch_localise_inv_subj.mâ.
- In the MATLAB command window, type âopen batch_localise_inv_subj_job.mâ. The cleared values have been replaced by â
<UNDEFINED>
â. We will provide these undefined values via the Matlab script. - In batch_localise_inv_subj.m replace nrun = X; with nrun = 1 and as inputs specify â{âPapMcbdspmeeg_run_01_sss.matâ}â and [-100 800].
- Finally, save and run the batch_localise_inv_subj.m script.
Review inverse solutions
- Review the inverse results by pressing the â3D Source Reconstructionâ button within the SPM window; re-load the âPapMcbdspmeeg_run_01_sss.matâ file.
- Press the âmipâ (âMaximum intensity projectionâ) button below the âInvertâ button. The top plot shows the evoked responses for the three conditions from the peak vertex, with the red line being the currently selected condition, here â1â for Famous faces (press the âconditionâ button to toggle through the other conditions). The bottom plot will show the maximal intensity projection at the time of the maximal activation.
- Take a look at the log model evidence and the percent variance explained.
- If you press âdisplayâ under the âWindowâ button, you can see a MIP for the time-frequency contrast limited to the 100-250ms, 10-20Hz window specified above; if you press the âdisplayâ under the âImageâ button, you will see a rendered version.
- Press the âpreviousâ button to select the previous inversion, which here corresponds to the MSP inversion. Press the âmipâ button again. You can change the time of interest or the vertex by entering new values into the window below the âmipâ button.
- Compare the model evidences and the percentage explained of the data variance of the MNM and the MSP solution. Which model explains the data better?
- Again, press âdisplayâ beneath the âWindowâ button and you should see results that are sparser and deeper inside the brain.
- You can further explore your inverse solutions by pressing the âRenderâ button of the â3D Source Reconstructionâ window. For example, you can play a movie of the source activity over time, look at the time courses of virtual electrodes (just first click on âvirtual electrodesâ and then on the area of the cortical surface you are interested in) or compare the observed and predicted sensor signals at different time points. You can toggle between models (MSP and MMN) and the conditions (1-Famous, 2-Unfamilar, 3 Scrambled) by pressing the buttons at the top of the render window.
We are at the end of the demo now. You can use the forward model under index â3â and play around with the inversion models, take a look at the induced activity or try to invert MEG and EEG data together.
SPM offers many other options for source reconstruction. Here we used two different distributed source reconstruction methods, but you could also use a Bayesian version of an Equivalent Current Dipole algorithm by selecting VB-ECD after clicking the âInvertâ button within the â3D Source Reconstructionâ window.
Literature
- http://www.fil.ion.ucl.ac.uk/spm/doc/manual.pdf
- Henson, R.N., Mattout, J., Phillips, C. and Friston, K.J. (2009). Selecting forward models for MEG source reconstruction using model-evidence. Neuroimage, 46, 168-176.
- Henson, R.N., Wakeman, D.G., Litvak, V. and Friston, K.J. (2011). A Parametric Empirical Bayesian framework for the EEG/MEG inverse problem: generative models for multisubject and multimodal integration. Frontiers in Human Neuroscience, 5, 76, 1-16. http://www.fil.ion.ucl.ac.uk/spm/doc/papers/HensonEtAl_FiHN_11_PEB_MEEG.pdf