Skip to content

Commit

Permalink
Tutorial: Add cMEM estimates in tutorials (#721)
Browse files Browse the repository at this point in the history
* Update `tutorial_epilepsy.m`
* Create `tutorial_BEst.m`
---------

Co-authored-by: rcassani <[email protected]>
  • Loading branch information
Edouard2laire and rcassani authored Nov 21, 2024
1 parent 9371938 commit a8784a4
Show file tree
Hide file tree
Showing 2 changed files with 241 additions and 0 deletions.
209 changes: 209 additions & 0 deletions toolbox/script/tutorial_BEst.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,209 @@
function tutorial_BEst(tutorial_dir, reports_dir)
% TUTORIAL_BEST: Script that runs all the Brain Entropy in space and time introduction tutorials.
%
% CORRESPONDING ONLINE TUTORIAL:
% https://neuroimage.usc.edu/brainstorm/Tutorials/TutBEst
%
% INPUTS:
% - tutorial_dir : Directory where the sample_introduction.zip file has been unzipped
% - reports_dir : Directory where to save the execution report (instead of displaying it)

% @=============================================================================
% This function is part of the Brainstorm software:
% https://neuroimage.usc.edu/brainstorm
%
% Copyright (c) University of Southern California & McGill University
% This software is distributed under the terms of the GNU General Public License
% as published by the Free Software Foundation. Further details on the GPLv3
% license can be found at http://www.gnu.org/copyleft/gpl.html.
%
% FOR RESEARCH PURPOSES ONLY. THE SOFTWARE IS PROVIDED "AS IS," AND THE
% UNIVERSITY OF SOUTHERN CALIFORNIA AND ITS COLLABORATORS DO NOT MAKE ANY
% WARRANTY, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO WARRANTIES OF
% MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE, NOR DO THEY ASSUME ANY
% LIABILITY OR RESPONSIBILITY FOR THE USE OF THIS SOFTWARE.
%
% For more information type "brainstorm license" at command prompt.
% =============================================================================@
%
% Author: Edouard Delaire, 2024
% Raymundo Cassani, 2024


% ===== FILES TO IMPORT =====
% Output folder for reports
if (nargin < 2) || isempty(reports_dir) || ~isdir(reports_dir)
reports_dir = [];
end
% Folder in which the Introduction tutorial dataset is unzipped (if needed)
if (nargin == 0) || isempty(tutorial_dir) || ~file_exist(tutorial_dir)
tutorial_dir = [];
end

% Re-inialize random number generator
if (bst_get('MatlabVersion') >= 712)
rng('default');
end


%% ===== VERIFY REQUIRED PROTOCOL =====
ProtocolName = 'TutorialIntroduction';
SubjectName = 'Subject01';

iProtocolIntroduction = bst_get('Protocol', ProtocolName);
if isempty(iProtocolIntroduction)
% Produce the Introduction protocol
tutorial_introduction(tutorial_dir, reports_dir)
else
% Select input protocol
gui_brainstorm('SetCurrentProtocol', iProtocolIntroduction);
end


%% ===== REQUIRED PLUGIN =====
% Install and Load Brain Entropy plugin
[isInstalled, errMsg] = bst_plugin('Install', 'brainentropy');
if ~isInstalled
error(errMsg);
end


%% ===== FIND FILES =====
% Process: Select recordings in: Subject01/S01_AEF_20131218_01_600Hz_notch
sFiles01 = bst_process('CallProcess', 'process_select_files_data', [], [], ...
'subjectname', SubjectName, ...
'condition', 'S01_AEF_20131218_01_600Hz_notch', ...
'includebad', 0);
% Process: Select file comments with tag: deviant
sFilesAvgDeviant01 = bst_process('CallProcess', 'process_select_tag', sFiles01, [], ...
'tag', 'Avg: deviant', ...
'search', 2, ... % Search the file comments
'select', 1); % Select only the files with the tag


%% ===== HEAD MODEL =====
disp([10 'BST> Head model' 10]);

% Process: Generate BEM surfaces
bst_process('CallProcess', 'process_generate_bem', [], [], ...
'subjectname', SubjectName, ...
'nscalp', 1922, ...
'nouter', 1922, ...
'ninner', 1922, ...
'thickness', 4, ...
'method', 'brainstorm');
% Process: Compute head model
bst_process('CallProcess', 'process_headmodel', sFilesAvgDeviant01, [], ...
'sourcespace', 1, ... % Cortex surface
'meg', 4, ... % OpenMEEG BEM
'openmeeg', struct(...
'BemSelect', [0, 0, 1], ...
'BemCond', [0.33, 0.0165, 0.33], ...
'BemNames', {{'Scalp', 'Skull', 'Brain'}}, ...
'BemFiles', {{}}, ...
'isAdjoint', 1, ...
'isAdaptative', 1, ...
'isSplit', 0, ...
'SplitLength', 4000));


%% ===== SOURCE ESTIMATION =====
% coherent Maximum Entropy on the Mean (cMEM)
disp([10 'BST> Source estimation using cMEM' 10]);

% Process: Compute sources: BEst
mem_option = be_pipelineoptions(be_main, 'cMEM');
mem_option.optional = struct_copy_fields(mem_option.optional, ...
struct(...
'TimeSegment', [0.05, 0.15], ...
'BaselineType', {{'within-data'}}, ...
'Baseline', [], ...
'BaselineHistory', {{'within'}}, ...
'BaselineSegment', [-0.1, 0], ...
'groupAnalysis', 0, ...
'display', 0));
sAvgSrcMEM = bst_process('CallProcess', 'process_inverse_mem', sFilesAvgDeviant01, [], ...
'comment', 'MEM', ...
'mem', struct('MEMpaneloptions', mem_option), ...
'sensortypes', 'MEG');
% Process: Snapshot: Sources (one time)
bst_process('CallProcess', 'process_snapshot', sAvgSrcMEM, [], ...
'target', 8, ... % Sources (one time)
'modality', 1, ... % MEG (All)
'orient', 1, ... % left
'time', 83.3*1e-3, ...
'threshold', 0, ...
'Comment', 'Average Deviant (cMEM)');


% wavelet Maximum Entropy on the Mean (wMEM)
disp([10 'BST> Source estimation using WMEM' 10]);

% Process: Compute sources: BEst
wMEM_options = be_pipelineoptions(be_main, 'wMEM');
wMEM_options.optional = struct_copy_fields(wMEM_options.optional, ...
struct(...
'TimeSegment', [0.05, 0.15], ...
'BaselineType', {{'within-data'}}, ...
'Baseline', [], ...
'BaselineHistory', {{'within'}}, ...
'BaselineSegment', [-0.1, 0], ...
'groupAnalysis', 0, ...
'display', 0));

% 1. Localizing only scale 4:
wMEM_options.wavelet.selected_scales = [4];
sAvgSrwMEM_scale4 = bst_process('CallProcess', 'process_inverse_mem', sFilesAvgDeviant01, [], ...
'comment', 'MEM', ...
'mem', struct( 'MEMpaneloptions', wMEM_options), ...
'sensortypes', 'MEG');
% Process: Snapshot: Sources (one time)
bst_process('CallProcess', 'process_snapshot', sAvgSrwMEM_scale4, [], ...
'target', 8, ... % Sources (one time)
'modality', 1, ... % MEG (All)
'orient', 1, ... % left
'time', 83.3*1e-3, ...
'threshold', 0, ...
'Comment', 'Average Deviant (wMEM - scale 4)');

% 2. Localizing only scale 5:
wMEM_options.wavelet.selected_scales = [5];
sAvgSrwMEM_scale5 = bst_process('CallProcess', 'process_inverse_mem', sFilesAvgDeviant01, [], ...
'comment', 'MEM', ...
'mem', struct( 'MEMpaneloptions', wMEM_options), ...
'sensortypes', 'MEG');
% Process: Snapshot: Sources (one time)
bst_process('CallProcess', 'process_snapshot', sAvgSrwMEM_scale5, [], ...
'target', 8, ... % Sources (one time)
'modality', 1, ... % MEG (All)
'orient', 1, ... % left
'time', 83.3*1e-3, ...
'threshold', 0, ...
'Comment', 'Average Deviant (wMEM - scale 5 )');

% 3. Localizing all scales:
wMEM_options.wavelet.selected_scales = [1,2,3,4,5];
sAvgSrwMEM_scaleAll = bst_process('CallProcess', 'process_inverse_mem', sFilesAvgDeviant01, [], ...
'comment', 'MEM', ...
'mem', struct( 'MEMpaneloptions', wMEM_options), ...
'sensortypes', 'MEG');
% Process: Snapshot: Sources (one time)
bst_process('CallProcess', 'process_snapshot', sAvgSrwMEM_scaleAll, [], ...
'target', 8, ... % Sources (one time)
'modality', 1, ... % MEG (All)
'orient', 1, ... % left
'time', 83.3*1e-3, ...
'threshold', 0, ...
'Comment', 'Average Deviant (wMEM - all scale)');


%% ===== SAVE REPORT =====
% Save and display report
ReportFile = bst_report('Save', []);
if ~isempty(reports_dir) && ~isempty(ReportFile)
bst_report('Export', ReportFile, reports_dir);
else
bst_report('Open', ReportFile);
end

disp([10 'BST> tutorial_BEst: Done.' 10]);
32 changes: 32 additions & 0 deletions toolbox/script/tutorial_epilepsy.m
Original file line number Diff line number Diff line change
Expand Up @@ -283,6 +283,38 @@ function tutorial_epilepsy(tutorial_dir, reports_dir)
'Comment', 'Average spike');


% Install and Load Brain Entropy plugin
[isInstalled, errMsg] = bst_plugin('Install', 'brainentropy');
if isInstalled
% Use default options of cMEM
mem_option = be_pipelineoptions(be_main, 'cMEM');
mem_option.optional = struct_copy_fields(mem_option.optional, ...
struct(...
'TimeSegment', [-0.30078, 0.5], ...
'BaselineType', {{'within-data'}}, ...
'Baseline', [], ...
'BaselineHistory', {{'within'}}, ...
'BaselineSegment', [-0.30078, -0.10156], ...
'groupAnalysis', 0, ...
'display', 0));
% Process: Compute sources: BEst
sAvgSrcMEM = bst_process('CallProcess', 'process_inverse_mem', sFilesAvg, [], ...
'comment', 'MEM', ...
'mem', struct('MEMpaneloptions', mem_option), ...
'sensortypes', 'EEG');
% Process: Snapshot: Sources (one time)
bst_process('CallProcess', 'process_snapshot', sAvgSrcMEM, [], ...
'target', 8, ... % Sources (one time)
'modality', 1, ... % MEG (All)
'orient', 3, ... % top
'time', 0, ...
'threshold', 30, ...
'Comment', 'Average spike (cMEM)');
else
bst_report('Error', [], [], errMsg);
end


% ===== SOURCE ANALYSIS: VOLUME =====
% Process: Compute head model
bst_process('CallProcess', 'process_headmodel', sFilesAvg, [], ...
Expand Down

0 comments on commit a8784a4

Please sign in to comment.