Resting state recordings from the OMEGA database

Authors: Francois Tadel, Guiomar Niso, Elizabeth Bock, Sylvain Baillet

This tutorial introduces how to download resting state recordings from the OMEGA database and process them into Brainstorm. The goal is to evaluate where in the brain is distributed the power in specific frequency bands at rest. This example includes only 5 healthy subjects, but the same pipeline can be applied over larger populations.

Note that the operations used here are not detailed, the goal of this tutorial is not to introduce Brainstorm to new users. For in-depth explanations of the interface and theoretical foundations, please refer to the introduction tutorials.

The processing pipeline is published in this article:
Niso G, Tadel F, Bock E, Cousineau M, Santos A, Baillet S, Brainstorm Pipeline Analysis of Resting-State Data from the Open MEG Archive, Frontiers in Neuroscience, Mar 2019

OMEGA database

The Open MEG Archive (OMEGA) is the fruit of a collaborative effort by the McConnell Brain Imaging Centre and the Université de Montréal to build a centralized repository in which to regroup MEG data for open dissemination. This continuously expanding repository also contains anatomical MRI volumes, demographic and questionnaire information, and eventually will feature other forms of electrophysiological data (e.g. EEG, field, and cell recordings). OMEGA features the technological framework for multi-centric data aggregation, and is amongst the largest freely available resting-state MEG datasets presently available. Website: https://www.mcgill.ca/bic/resources/omega

You are free to use all data in OMEGA for research purposes; however, we ask that you please cite the following reference in your publications if you have used data from OMEGA:

Niso G, Rogers C, Moreau JT, Chen LY, Madjar C, Das S, Bock E, Tadel F, Evans A, Jolicoeur P, Baillet S
OMEGA: The Open MEG Archive, Neuroimage, 2015.

Presentation of the experiment

Experiment

MEG acquisition

Head shape and fiducial points

3D digitization using a Polhemus Fastrak device driven by Brainstorm. The .pos files contain:

Subject anatomy

Note that while the subjects shared here have their anatomy already processed with FreeSurfer to simplify the tutorial, this is not currently the case for data shared in OMEGA. You would therefore need to generate surfaces yourself as explained in the anatomy tutorial.

Download and installation

BIDS specification

The Brain Imaging Data Structure (BIDS) standard for neuroimaging data organization was first established for MRI and fMRI (Gorgolewski, 2016). It is based on simple file formats and folder structures that can readily expand to additional data modalities.

An extension of the BIDS for MEG datasets has now been published. The data distributed with this tutorial is organized following this prototype (which may change in the future). One objective is that analysis pipelines designed with major analysis tools (such as Brainstorm, FieldTrip, MNE, SPM and others) can be readily applied without requiring software or pipeline redevelopments.

Brainstorm includes a set of functions that can load automatically any dataset following the BIDS specification. The files that are read by Brainstorm are the following:

OMEGA_RestingState_sample/

Import the dataset

Refine the MEG-MRI registration

Pre-processing

Artifact cleaning

We will now run the automatic procedure for cleaning the heartbeats, as described in the introduction tutorials (detection, SSP). The results we want to illustrate here are robust enough, the recordings do not need to be processed any further. If you want to improve the quality of the data with more manual cleaning (blinks, saccades, movements, bad segments), please refer to the introduction tutorials.

Source estimation

Power maps

We have now access to continuous recordings in source space, we can estimate the average power over this resting period, for various frequency bands. We propose here to compute the PSD from the first 100s instead of the full 600s rest recordings, it is faster and leads to very similar results. For more stable estimates, you can compute the the PSD over all the data available instead.

Scripting

The following script from the Brainstorm distribution reproduces the analysis presented in this tutorial page: brainstorm3/toolbox/script/tutorial_omega.m

1 function tutorial_omega(BidsDir) 2 % TUTORIAL_OMEGA: Script that reproduces the results of the online tutorial "Resting state and OMEGA database". 3 % 4 % CORRESPONDING ONLINE TUTORIALS: 5 % https://neuroimage.usc.edu/brainstorm/Tutorials/RestingOmega 6 % 7 % INPUTS: 8 % tutorial_dir: Directory where the OMEGA_RestingState_sample.tar file has been unzipped 9 10 % @============================================================================= 11 % This function is part of the Brainstorm software: 12 % https://neuroimage.usc.edu/brainstorm 13 % 14 % Copyright (c) University of Southern California & McGill University 15 % This software is distributed under the terms of the GNU General Public License 16 % as published by the Free Software Foundation. Further details on the GPLv3 17 % license can be found at http://www.gnu.org/copyleft/gpl.html. 18 % 19 % FOR RESEARCH PURPOSES ONLY. THE SOFTWARE IS PROVIDED "AS IS," AND THE 20 % UNIVERSITY OF SOUTHERN CALIFORNIA AND ITS COLLABORATORS DO NOT MAKE ANY 21 % WARRANTY, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO WARRANTIES OF 22 % MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE, NOR DO THEY ASSUME ANY 23 % LIABILITY OR RESPONSIBILITY FOR THE USE OF THIS SOFTWARE. 24 % 25 % For more information type "brainstorm license" at command prompt. 26 % =============================================================================@ 27 % 28 % Author: Francois Tadel, 2016-2018 29 30 31 %% ===== FILES TO IMPORT ===== 32 % You have to specify the folder in which the tutorial dataset is unzipped 33 if (nargin == 0) || isempty(BidsDir) || ~file_exist(BidsDir) 34 error('The first argument must be the full path to the tutorial dataset folder.'); 35 end 36 37 38 %% ===== CREATE PROTOCOL ===== 39 % The protocol name has to be a valid folder name (no spaces, no weird characters...) 40 ProtocolName = 'TutorialOmega'; 41 % Start brainstorm without the GUI 42 if ~brainstorm('status') 43 brainstorm nogui 44 end 45 % Delete existing protocol 46 gui_brainstorm('DeleteProtocol', ProtocolName); 47 % Create new protocol 48 gui_brainstorm('CreateProtocol', ProtocolName, 0, 0); 49 % Start a new report 50 bst_report('Start'); 51 52 53 %% ===== IMPORT BIDS DATASET ===== 54 % Process: Import BIDS dataset 55 sFilesRaw = bst_process('CallProcess', 'process_import_bids', [], [], ... 56 'bidsdir', {BidsDir, 'BIDS'}, ... 57 'nvertices', 15000, ... 58 'channelalign', 0); 59 60 % Process: Remove head points 61 sFilesRaw = bst_process('CallProcess', 'process_headpoints_remove', sFilesRaw, [], ... 62 'zlimit', 0); 63 64 % Process: Refine registration 65 sFilesRaw = bst_process('CallProcess', 'process_headpoints_refine', sFilesRaw, []); 66 67 % Process: Convert to continuous (CTF): Continuous 68 sFilesRaw = bst_process('CallProcess', 'process_ctf_convert', sFilesRaw, [], ... 69 'rectype', 2); % Continuous 70 71 72 %% ===== PRE-PROCESSING ===== 73 % % Process: Power spectrum density (Welch) 74 % sFilesPsdBefore = bst_process('CallProcess', 'process_psd', sFilesRaw, [], ... 75 % 'timewindow', [], ... 76 % 'win_length', 4, ... 77 % 'win_overlap', 50, ... 78 % 'sensortypes', 'MEG, EEG', ... 79 % 'edit', struct(... 80 % 'Comment', 'Power', ... 81 % 'TimeBands', [], ... 82 % 'Freqs', [], ... 83 % 'ClusterFuncTime', 'none', ... 84 % 'Measure', 'power', ... 85 % 'Output', 'all', ... 86 % 'SaveKernel', 0)); 87 88 % Process: Notch filter: 60Hz 120Hz 180Hz 240Hz 300Hz 89 sFilesNotch = bst_process('CallProcess', 'process_notch', sFilesRaw, [], ... 90 'freqlist', [60, 120, 180, 240, 300], ... 91 'sensortypes', 'MEG, EEG', ... 92 'read_all', 1); 93 94 % Process: High-pass:0.3Hz 95 sFilesBand = bst_process('CallProcess', 'process_bandpass', sFilesNotch, [], ... 96 'sensortypes', 'MEG, EEG', ... 97 'highpass', 0.3, ... 98 'lowpass', 0, ... 99 'attenuation', 'strict', ... % 60dB 100 'mirror', 0, ... 101 'useold', 0, ... 102 'read_all', 1); 103 104 % Process: Power spectrum density (Welch) 105 sFilesPsdAfter = bst_process('CallProcess', 'process_psd', sFilesBand, [], ... 106 'timewindow', [0 100], ... 107 'win_length', 4, ... 108 'win_overlap', 50, ... 109 'sensortypes', 'MEG, EEG', ... 110 'edit', struct(... 111 'Comment', 'Power', ... 112 'TimeBands', [], ... 113 'Freqs', [], ... 114 'ClusterFuncTime', 'none', ... 115 'Measure', 'power', ... 116 'Output', 'all', ... 117 'SaveKernel', 0)); 118 119 % Process: Snapshot: Frequency spectrum 120 bst_process('CallProcess', 'process_snapshot', sFilesPsdAfter, [], ... 121 'target', 10, ... % Frequency spectrum 122 'modality', 1); % MEG (All) 123 124 % Process: Delete folders 125 bst_process('CallProcess', 'process_delete', [sFilesRaw, sFilesNotch], [], ... 126 'target', 2); % Delete folders 127 128 129 %% ===== ARTIFACT CLEANING ===== 130 % Process: Select data files in: */* 131 sFilesBand = bst_process('CallProcess', 'process_select_files_data', [], [], ... 132 'subjectname', 'All'); 133 134 % Process: Select file names with tag: task-rest 135 sFilesRest = bst_process('CallProcess', 'process_select_tag', sFilesBand, [], ... 136 'tag', 'task-rest', ... 137 'search', 1, ... % Search the file names 138 'select', 1); % Select only the files with the tag 139 140 % Process: Detect heartbeats 141 bst_process('CallProcess', 'process_evt_detect_ecg', sFilesRest, [], ... 142 'channelname', 'ECG', ... 143 'timewindow', [], ... 144 'eventname', 'cardiac'); 145 146 % Process: SSP ECG: cardiac 147 bst_process('CallProcess', 'process_ssp_ecg', sFilesRest, [], ... 148 'eventname', 'cardiac', ... 149 'sensortypes', 'MEG', ... 150 'usessp', 1, ... 151 'select', 1); 152 153 % Process: Snapshot: Sensors/MRI registration 154 bst_process('CallProcess', 'process_snapshot', sFilesRest, [], ... 155 'target', 1, ... % Sensors/MRI registration 156 'modality', 1, ... % MEG (All) 157 'orient', 1); % left 158 159 % Process: Snapshot: SSP projectors 160 bst_process('CallProcess', 'process_snapshot', sFilesRest, [], ... 161 'target', 2, ... % SSP projectors 162 'modality', 1); % MEG (All) 163 164 165 %% ===== SOURCE ESTIMATION ===== 166 % Process: Select file names with tag: task-rest 167 sFilesNoise = bst_process('CallProcess', 'process_select_tag', sFilesBand, [], ... 168 'tag', 'task-noise', ... 169 'search', 1, ... % Search the file names 170 'select', 1); % Select only the files with the tag 171 172 % Process: Compute covariance (noise or data) 173 bst_process('CallProcess', 'process_noisecov', sFilesNoise, [], ... 174 'baseline', [], ... 175 'sensortypes', 'MEG', ... 176 'target', 1, ... % Noise covariance (covariance over baseline time window) 177 'dcoffset', 1, ... % Block by block, to avoid effects of slow shifts in data 178 'identity', 0, ... 179 'copycond', 1, ... 180 'copysubj', 1, ... 181 'copymatch', 1, ... 182 'replacefile', 1); % Replace 183 184 % Process: Compute head model 185 bst_process('CallProcess', 'process_headmodel', sFilesRest, [], ... 186 'sourcespace', 1, ... % Cortex surface 187 'meg', 3); % Overlapping spheres 188 189 % Process: Compute sources [2018] 190 sSrcRest = bst_process('CallProcess', 'process_inverse_2018', sFilesRest, [], ... 191 'output', 2, ... % Kernel only: one per file 192 'inverse', struct(... 193 'Comment', 'dSPM: MEG', ... 194 'InverseMethod', 'minnorm', ... 195 'InverseMeasure', 'dspm2018', ... 196 'SourceOrient', {{'fixed'}}, ... 197 'Loose', 0.2, ... 198 'UseDepth', 1, ... 199 'WeightExp', 0.5, ... 200 'WeightLimit', 10, ... 201 'NoiseMethod', 'reg', ... 202 'NoiseReg', 0.1, ... 203 'SnrMethod', 'fixed', ... 204 'SnrRms', 1e-06, ... 205 'SnrFixed', 3, ... 206 'ComputeKernel', 1, ... 207 'DataTypes', {{'MEG'}})); 208 209 210 %% ===== POWER MAPS ===== 211 % Process: Power spectrum density (Welch) 212 sSrcPsd = bst_process('CallProcess', 'process_psd', sSrcRest, [], ... 213 'timewindow', [0, 100], ... 214 'win_length', 4, ... 215 'win_overlap', 50, ... 216 'clusters', {}, ... 217 'scoutfunc', 1, ... % Mean 218 'edit', struct(... 219 'Comment', 'Power,FreqBands', ... 220 'TimeBands', [], ... 221 'Freqs', {{'delta', '2, 4', 'mean'; 'theta', '5, 7', 'mean'; 'alpha', '8, 12', 'mean'; 'beta', '15, 29', 'mean'; 'gamma1', '30, 59', 'mean'; 'gamma2', '60, 90', 'mean'}}, ... 222 'ClusterFuncTime', 'none', ... 223 'Measure', 'power', ... 224 'Output', 'all', ... 225 'SaveKernel', 0)); 226 227 % Process: Spectrum normalization 228 sSrcPsdNorm = bst_process('CallProcess', 'process_tf_norm', sSrcPsd, [], ... 229 'normalize', 'relative', ... % Relative power (divide by total power) 230 'overwrite', 0); 231 232 % Process: Project on default anatomy: surface 233 sSrcPsdProj = bst_process('CallProcess', 'process_project_sources', sSrcPsdNorm, [], ... 234 'headmodeltype', 'surface'); % Cortex surface 235 236 % Process: Spatial smoothing (3.00) 237 sSrcPsdProj = bst_process('CallProcess', 'process_ssmooth_surfstat', sSrcPsdProj, [], ... 238 'fwhm', 3, ... 239 'overwrite', 1); 240 241 % Process: Average: Everything 242 sSrcPsdAvg = bst_process('CallProcess', 'process_average', sSrcPsdProj, [], ... 243 'avgtype', 1, ... % Everything 244 'avg_func', 1, ... % Arithmetic average: mean(x) 245 'weighted', 0, ... 246 'matchrows', 0, ... 247 'iszerobad', 0); 248 249 % Screen capture of final result 250 hFig = view_surface_data([], sSrcPsdAvg.FileName); 251 set(hFig, 'Position', [200 200 200 200]); 252 hFigContact = view_contactsheet(hFig, 'freq', 'fig'); 253 bst_report('Snapshot', hFigContact, sSrcPsdAvg.FileName, 'Power'); 254 close([hFig, hFigContact]); 255 256 % Save and display report 257 ReportFile = bst_report('Save', []); 258 bst_report('Open', ReportFile); 259 260 261 262

Additional documentation

Articles

Forum





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


Tutorials/RestingOmega (last edited 2022-11-28 10:28:33 by FrancoisTadel)