MEG resting state and phase-amplitude coupling

Authors: Thomas Donoghue, Soheila Samiee, Elizabeth Bock, Esther Florin, Francois Tadel, Sylvain Baillet

This tutorial features a few pre-processing steps and the calculation of phase-amplitude coupling measures on resting state MEG recordings. It is based on a eyes open resting recordings of one subject recorded at the Montreal Neurological Institute in 2012 with a CTF MEG 275 system.

License

This tutorial dataset (MEG and MRI data) remains a property of the MEG Lab, McConnell Brain Imaging Center, Montreal Neurological Institute, McGill University, Canada. Its use and transfer outside the Brainstorm tutorial, e.g. for research purposes, is prohibited without written consent from the MEG Lab.

If you reference this dataset in your publications, please acknowledge its authors (Elizabeth Bock, Esther Florin, Francois Tadel and Sylvain Baillet) and cite Brainstorm as indicated on the website. For questions, please contact us through the forum.

Presentation of the experiment

Experiment

  • One subject
  • Two runs of 10 min of resting state recordings
  • Eyes open

MEG acquisition

  • Acquisition at 2400Hz, with a CTF 275 system, subject in seating position

  • Recorded at the Montreal Neurological Institute in 2012
  • Anti-aliasing low-pass filter at 600Hz, files saved with the 3rd order gradient
  • Recorded channels (301):
    • 26 MEG reference sensors (#2-#27)
    • 272 MEG axial gradiometers (#28-#299)
    • 1 ECG bipolar (#300)
    • 2 vertical EOG bipolar (#301)
    • 1 Unused channels (#1)
  • 3 datasets:
    • subj002_spontaneous_20111102_01_AUX.ds: Run #1, 10min resting state

    • subj002_spontaneous_20111102_02_AUX.ds: Run #2, 10min resting state

    • subj002_noise_20111104_02.ds: Empty room recordings, 100s

Head shape and fiducial points

  • 3D digitization using a Polhemus Fastrak device (subj002_11022011.pos)

  • Position of the anatomical references used for these recordings:
    The position of the CTF head localization coils, pictures available in file fiducials.jpg

  • Around 100 head points distributed on the hard parts of the head (no soft tissues)

Subject anatomy

  • Subject with 1.5T MRI
  • Processed with FreeSurfer 5.2

Download and installation

Import the anatomy

Access the recordings

Resting state recordings: 10min

Empty room recordings: 90s

Pre-processing

All data should be pre-processed and checked for artifacts prior to doing analyses such as PAC (including marking bad segments and correcting for artifacts). For the purposes of this tutorial, we will correct for blinks and heartbeats with SSPs but will not go through marking out bad sections. When using your own data reviewing the raw data for bad segments and using clean data is of the utmost importance.

Bad segments

Some noisy segments have been marked manually. We are going to import these segments now:

Signal Space Projection (SSP) is a method for projecting the recordings away from stereotyped artifacts, such as eye blinks and heartbeats.

Power line contamination

Source estimation

We need now to calculate a source model for the resting state recordings, using a noise covariance matrix calculated from the noise recordings.

Head model

Noise covariance

Inverse model

Scouts

Phase-amplitude coupling

We are now ready to run the PAC analysis on the source signals. This PAC function in Brainstorm is not time resolved, but will analyze the given time series for any stable occurrence of PAC over a time segment you give it. For more information about the PAC measure used here, please refer to the online tutorial
Phase-amplitude coupling.

PAC estimation

Visual exploration of the comodulogram

Scout functions

Canolty maps

Introduction to the method

Canolty maps are a type of time-frequency decomposition that offer another way to visualize the data and serve as a complimentary tool to visualize and assess phase-amplitude coupling. The process lines up the data to a specific low frequency so as to visualize what happens in the power spectrum related to the phase of the low frequency. Currently there are no significance tests within Brainstorm that can give a measure if PAC is significant in a given time series, but the Canolty maps provide an important way to verify and corroborate the results of the PAC process. We name these maps after the author of the paper in which they were first introduced:

Canolty RT, Edwards E, Dalal SS, Soltani M, Nagarajan SS, Kirsch HE, Berger MS, Barbaro NM, Knight RT
High gamma power is phase-locked to theta oscillations in human neocortex
Science, 2006 Sep 15;313(5793):1626-8.

The procedure to obtain a Canolty map for a signal is the following:

Computation

Scout functions

Other frequencies

Two inputs

Scripting

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

1 function tutorial_resting(tutorial_dir) 2 % TUTORIAL_RESTING: Script that reproduces the results of the online tutorial "Resting state MEG". 3 % 4 % DESCRIPTION: 5 % This tutorial is based on two blocks of 10 minutes of resting state eyes open recorded 6 % at the Montreal Neurological Institute in 2011 with a CTF MEG 275 system. 7 % 8 % INPUTS: 9 % tutorial_dir: Directory where the sample_resting.zip file has been unzipped 10 11 % @============================================================================= 12 % This function is part of the Brainstorm software: 13 % https://neuroimage.usc.edu/brainstorm 14 % 15 % Copyright (c) University of Southern California & McGill University 16 % This software is distributed under the terms of the GNU General Public License 17 % as published by the Free Software Foundation. Further details on the GPLv3 18 % license can be found at http://www.gnu.org/copyleft/gpl.html. 19 % 20 % FOR RESEARCH PURPOSES ONLY. THE SOFTWARE IS PROVIDED "AS IS," AND THE 21 % UNIVERSITY OF SOUTHERN CALIFORNIA AND ITS COLLABORATORS DO NOT MAKE ANY 22 % WARRANTY, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO WARRANTIES OF 23 % MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE, NOR DO THEY ASSUME ANY 24 % LIABILITY OR RESPONSIBILITY FOR THE USE OF THIS SOFTWARE. 25 % 26 % For more information type "brainstorm license" at command prompt. 27 % =============================================================================@ 28 % 29 % Author: Francois Tadel, 2013-2016 30 31 32 % ===== FILES TO IMPORT ===== 33 % You have to specify the folder in which the tutorial dataset is unzipped 34 if (nargin == 0) || isempty(tutorial_dir) || ~file_exist(tutorial_dir) 35 error('The first argument must be the full path to the tutorial dataset folder.'); 36 end 37 % Build the files name of the files to import 38 AnatDir = fullfile(tutorial_dir, 'sample_resting', 'Anatomy'); 39 Run1File = fullfile(tutorial_dir, 'sample_resting', 'Data', 'subj002_spontaneous_20111102_01_AUX.ds'); 40 Run2File = fullfile(tutorial_dir, 'sample_resting', 'Data', 'subj002_spontaneous_20111102_02_AUX.ds'); 41 NoiseFile = fullfile(tutorial_dir, 'sample_resting', 'Data', 'subj002_noise_20111104_02.ds'); 42 Event1File = fullfile(Run1File, 'events_BAD.mat'); 43 Event2File = fullfile(Run2File, 'events_BAD.mat'); 44 % Check if the folder contains the required files 45 if ~file_exist(Run1File) 46 error(['The folder ' tutorial_dir ' does not contain the folder from the file sample_resting.zip.']); 47 end 48 49 % ===== CREATE PROTOCOL ===== 50 % The protocol name has to be a valid folder name (no spaces, no weird characters...) 51 ProtocolName = 'TutorialResting'; 52 % Start brainstorm without the GUI 53 if ~brainstorm('status') 54 brainstorm nogui 55 end 56 % Delete existing protocol 57 gui_brainstorm('DeleteProtocol', ProtocolName); 58 % Create new protocol 59 gui_brainstorm('CreateProtocol', ProtocolName, 0, 0); 60 % Start a new report 61 bst_report('Start'); 62 63 64 % ===== ANATOMY ===== 65 % Subject name 66 SubjectName = 'Subject02'; 67 % Process: Import FreeSurfer folder 68 bst_process('CallProcess', 'process_import_anatomy', [], [], ... 69 'subjectname', SubjectName, ... 70 'mrifile', {AnatDir, 'FreeSurfer'}, ... 71 'nvertices', 15000, ... 72 'nas', [128, 225, 135], ... 73 'lpa', [ 54, 115, 107], ... 74 'rpa', [204, 115, 99]); 75 76 77 % ===== LINK CONTINUOUS FILE ===== 78 % Process: Create link to raw file (Run1) 79 sFileRun1 = bst_process('CallProcess', 'process_import_data_raw', [], [], ... 80 'subjectname', SubjectName, ... 81 'datafile', {Run1File, 'CTF'}, ... 82 'channelalign', 1); 83 84 % Process: Events: Import from file (Run1) 85 bst_process('CallProcess', 'process_evt_import', sFileRun1, [], ... 86 'evtfile', {Event1File, 'BST'}); 87 88 % Link noise recordings (Noise) 89 sFileNoise = bst_process('CallProcess', 'process_import_data_raw', [], [], ... 90 'subjectname', SubjectName, ... 91 'datafile', {NoiseFile, 'CTF'}, ... 92 'channelalign', 0); 93 94 % Process: Snapshot: Sensors/MRI registration (Run1) 95 bst_process('CallProcess', 'process_snapshot', sFileRun1, [], ... 96 'target', 1, ... % Sensors/MRI registration 97 'modality', 1, ... % MEG (All) 98 'orient', 1, ... % left 99 'Comment', 'MEG/MRI Registration'); 100 101 % Process: Create link to raw file 102 % sFilesRun2 = bst_process('CallProcess', 'process_import_data_raw', [], [], ... 103 % 'subjectname', SubjectName, ... 104 % 'datafile', {Run2File, 'CTF'}, ... 105 % 'channelalign', 1); 106 % Process: Events: Import from file 107 % bst_process('CallProcess', 'process_evt_import', sFilesRun2, [], ... 108 % 'evtfile', {Event2File, 'BST'}); 109 110 111 % ===== CORRECT BLINKS AND HEARTBEATS ===== 112 % Process: Detect heartbeats (Run1) 113 sFileRun1 = bst_process('CallProcess', 'process_evt_detect_ecg', sFileRun1, [], ... 114 'channelname', 'EEG057', ... 115 'timewindow', [], ... 116 'eventname', 'cardiac'); 117 118 % Process: Detect eye blinks (Run1) 119 sFileRun1 = bst_process('CallProcess', 'process_evt_detect_eog', sFileRun1, [], ... 120 'channelname', 'EEG058', ... 121 'timewindow', [], ... 122 'eventname', 'blink'); 123 124 % Process: Remove simultaneous (Run1) 125 sFileRun1 = bst_process('CallProcess', 'process_evt_remove_simult', sFileRun1, [], ... 126 'remove', 'cardiac', ... 127 'target', 'blink', ... 128 'dt', 0.25, ... 129 'rename', 0); 130 131 % Process: SSP ECG: cardiac (Run1) 132 sFileRun1 = bst_process('CallProcess', 'process_ssp_ecg', sFileRun1, [], ... 133 'eventname', 'cardiac', ... 134 'sensortypes', 'MEG', ... 135 'usessp', 1, ... 136 'select', 1); 137 138 % Process: SSP EOG: blink (Run1) 139 sFileRun1 = bst_process('CallProcess', 'process_ssp_eog', sFileRun1, [], ... 140 'eventname', 'blink', ... 141 'sensortypes', 'MEG', ... 142 'usessp', 1, ... 143 'select', 1); 144 145 % Process: Snapshot: SSP projectors (Run1) 146 bst_process('CallProcess', 'process_snapshot', sFileRun1, [], ... 147 'target', 2, ... % SSP projectors 148 'Comment', 'SSP projectors'); 149 150 151 % ===== SOURCE MODELING ===== 152 % Process: Compute head model (Run1) 153 bst_process('CallProcess', 'process_headmodel', sFileRun1, [], ... 154 'comment', '', ... 155 'sourcespace', 1, ... 156 'meg', 3); % Overlapping spheres 157 158 % Process: Compute noise covariance (Noise) 159 bst_process('CallProcess', 'process_noisecov', sFileNoise, [], ... 160 'baseline', [], ... 161 'sensortypes', 'MEG', ... % Noise covariance (covariance over baseline time window) 162 'dcoffset', 1, ... 163 'identity', 0, ... 164 'copycond', 1, ... 165 'copysubj', 0, ... 166 'replacefile', 1); % Replace 167 168 % Process: Snapshot: Noise covariance 169 bst_process('CallProcess', 'process_snapshot', sFileNoise, [], ... 170 'target', 3, ... % Noise covariance 171 'Comment', 'Noise covariance'); 172 173 % Process: Compute sources (Run1) 174 sFileRun1Source = bst_process('CallProcess', 'process_inverse', sFileRun1, [], ... 175 'comment', '', ... 176 'method', 1, ... % Minimum norm estimates (wMNE) 177 'wmne', struct(... 178 'NoiseCov', [], ... 179 'InverseMethod', 'wmne', ... 180 'ChannelTypes', {{}}, ... 181 'SNR', 3, ... 182 'diagnoise', 0, ... 183 'SourceOrient', {{'fixed'}}, ... 184 'loose', 0.2, ... 185 'depth', 1, ... 186 'weightexp', 0.5, ... 187 'weightlimit', 10, ... 188 'regnoise', 1, ... 189 'magreg', 0.1, ... 190 'gradreg', 0.1, ... 191 'eegreg', 0.1, ... 192 'ecogreg', 0.1, ... 193 'seegreg', 0.1, ... 194 'fMRI', [], ... 195 'fMRIthresh', [], ... 196 'fMRIoff', 0.1, ... 197 'pca', 1), ... 198 'sensortypes', 'MEG', ... 199 'output', 2); % Kernel only: one per file 200 201 % ===== RESTING STATE PIPELINE ===== 202 % Process: Phase-amplitude coupling 203 sFilePac = bst_process('CallProcess', 'process_pac', sFileRun1Source, [], ... 204 'timewindow', [400, 600], ... 205 'scouts', {'Brodmann', {'V1 R'}}, ... 206 'scoutfunc', 1, ... % Mean 207 'scouttime', 1, ... % Before 208 'nesting', [2, 14], ... 209 'nested', [40, 150], ... 210 'numfreqs', 0, ... 211 'parallel', 0, ... 212 'ismex', 1, ... 213 'max_block_size', 1, ... 214 'avgoutput', 0, ... 215 'savemax', 0); 216 217 % Process: Canolty maps 218 sFileCanolty = bst_process('CallProcess', 'process_canoltymap', sFileRun1Source, [], ... 219 'timewindow', [400, 600], ... 220 'scouts', {'Brodmann', {'V1 R'}}, ... 221 'scoutfunc', 5, ... % All 222 'scouttime', 2, ... 223 'epochtime', [-0.5, 0.5], ... 224 'lowfreq', 11, ... 225 'max_block_size', 100, ... 226 'save_erp', 1); 227 228 % Process: Canolty maps (FileB=MaxPAC) 229 sFileCanolty2 = bst_process('CallProcess', 'process_canoltymap2', sFileRun1Source, sFilePac, ... 230 'timewindow', [400, 600], ... 231 'scouts', {'Brodmann', {'V1 R'}}, ... 232 'scoutfunc', 5, ... % All 233 'scouttime', 1, ... 234 'epochtime', [-0.5, 0.5], ... 235 'max_block_size', 100, ... 236 'save_erp', 1); 237 238 % Display results 239 hFig1 = view_pac(sFilePac.FileName); 240 hFig2 = view_timefreq(sFileCanolty.FileName); 241 hFig3 = view_timefreq(sFileCanolty2.FileName); 242 pause(1); 243 close([hFig1 hFig2 hFig3]); 244 245 % Save and display report 246 ReportFile = bst_report('Save', sFileCanolty2); 247 bst_report('Open', ReportFile); 248 249 250 251 252





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


Tutorials/Resting (last edited 2022-06-24 10:20:08 by FrancoisTadel)