Created
May 26, 2023 20:10
-
-
Save rcassani/3ebbcebdd7e12aa11b1054a9d91f846e to your computer and use it in GitHub Desktop.
Connectivity tests with CMC tutorial data
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| % Script to run connectivity (magnitude square coherence) on data from the | |
| % Corticomuscular coherence tutorial: | |
| % https://neuroimage.usc.edu/brainstorm/Tutorials/CorticomuscularCoherence | |
| % | |
| % This is done in order to verify the difference approaches for PCA | |
| % provided in the ScoutPca branch (by Marc) | |
| %% Parameters | |
| % Brainstorm version and Protocol | |
| ProtocolName = 'TutorialCMC'; | |
| SubjectName = 'Subject01'; | |
| ConditionName = 'SubjectCMC_test_pca_8afbc14'; | |
| GitBranchName = 'ScoutPca'; | |
| BrainstormPath = '/export02/data/rcassani/Documents/mcgill/brainstorm/brainstorm3-dev'; | |
| % Coherence | |
| emgChannel = 'EMGlft'; % Name of EMG channel | |
| atlasScouts = {'Schaefer_100_17net', {'Background+FreeSurfer_Defined_Medial_Wall L', 'Background+FreeSurfer_Defined_Medial_Wall R', 'ContA_IPS_1 L', 'ContA_IPS_1 R', 'ContA_PFCl_1 L', 'ContA_PFCl_1 R', 'ContA_PFCl_2 L', 'ContA_PFCl_2 R', 'ContB_IPL_1 R', 'ContB_PFCld_1 R', 'ContB_PFClv_1 L', 'ContB_PFClv_1 R', 'ContB_Temp_1 R', 'ContC_Cingp_1 L', 'ContC_Cingp_1 R', 'ContC_pCun_1 L', 'ContC_pCun_1 R', 'ContC_pCun_2 L', 'DefaultA_IPL_1 R', 'DefaultA_PFCd_1 L', 'DefaultA_PFCd_1 R', 'DefaultA_PFCm_1 L', 'DefaultA_PFCm_1 R', 'DefaultA_pCunPCC_1 L', 'DefaultA_pCunPCC_1 R', 'DefaultB_IPL_1 L', 'DefaultB_PFCd_1 L', 'DefaultB_PFCd_1 R', 'DefaultB_PFCl_1 L', 'DefaultB_PFCv_1 L', 'DefaultB_PFCv_1 R', 'DefaultB_PFCv_2 L', 'DefaultB_PFCv_2 R', 'DefaultB_Temp_1 L', 'DefaultB_Temp_2 L', 'DefaultC_PHC_1 L', 'DefaultC_PHC_1 R', 'DefaultC_Rsp_1 L', 'DefaultC_Rsp_1 R', 'DorsAttnA_ParOcc_1 L', 'DorsAttnA_ParOcc_1 R', 'DorsAttnA_SPL_1 L', 'DorsAttnA_SPL_1 R', 'DorsAttnA_TempOcc_1 L', 'DorsAttnA_TempOcc_1 R', 'DorsAttnB_FEF_1 L', 'DorsAttnB_FEF_1 R', 'DorsAttnB_PostC_1 L', 'DorsAttnB_PostC_1 R', 'DorsAttnB_PostC_2 L', 'DorsAttnB_PostC_2 R', 'DorsAttnB_PostC_3 L', 'LimbicA_TempPole_1 L', 'LimbicA_TempPole_1 R', 'LimbicA_TempPole_2 L', 'LimbicB_OFC_1 L', 'LimbicB_OFC_1 R', 'SalVentAttnA_FrMed_1 L', 'SalVentAttnA_FrMed_1 R', 'SalVentAttnA_Ins_1 L', 'SalVentAttnA_Ins_1 R', 'SalVentAttnA_Ins_2 L', 'SalVentAttnA_ParMed_1 L', 'SalVentAttnA_ParMed_1 R', 'SalVentAttnA_ParOper_1 L', 'SalVentAttnA_ParOper_1 R', 'SalVentAttnB_IPL_1 R', 'SalVentAttnB_PFCl_1 L', 'SalVentAttnB_PFCl_1 R', 'SalVentAttnB_PFCmp_1 L', 'SalVentAttnB_PFCmp_1 R', 'SomMotA_1 L', 'SomMotA_1 R', 'SomMotA_2 L', 'SomMotA_2 R', 'SomMotA_3 R', 'SomMotA_4 R', 'SomMotB_Aud_1 L', 'SomMotB_Aud_1 R', 'SomMotB_Cent_1 L', 'SomMotB_Cent_1 R', 'SomMotB_S2_1 L', 'SomMotB_S2_1 R', 'SomMotB_S2_2 L', 'SomMotB_S2_2 R', 'TempPar_1 L', 'TempPar_1 R', 'TempPar_2 R', 'TempPar_3 R', 'VisCent_ExStr_1 L', 'VisCent_ExStr_1 R', 'VisCent_ExStr_2 L', 'VisCent_ExStr_2 R', 'VisCent_ExStr_3 L', 'VisCent_ExStr_3 R', 'VisCent_Striate_1 L', 'VisPeri_ExStrInf_1 L', 'VisPeri_ExStrInf_1 R', 'VisPeri_ExStrSup_1 L', 'VisPeri_ExStrSup_1 R', 'VisPeri_StriCal_1 L', 'VisPeri_StriCal_1 R'}}; | |
| cohMeasure = 'mscohere'; % Magnitude-squared Coherence|C|^2 = |Cxy|^2/(Cxx*Cyy) | |
| winLength = 0.5; % 500ms | |
| overlap = 50; % 50% | |
| maxFreq = 80; % 80Hz | |
| %% Start | |
| startTime = datestr(datetime('now','Format','d-MMM-y HH:mm:ss')); | |
| % Verify that we are in the proper Brainstorm branch | |
| cd(BrainstormPath); | |
| [status,cmdout] = system('git branch --show-current'); | |
| if ~strcmpi(GitBranchName, cmdout(1:end-1)) | |
| error('Incorrect Git branch') | |
| end | |
| % Start Brainstorm | |
| brainstorm; | |
| % Select protocol | |
| iProtocol = bst_get('Protocol', ProtocolName); | |
| if isempty(iProtocol) | |
| error(['Unknown protocol: ' ProtocolName]); | |
| end | |
| if (iProtocol ~= bst_get('iProtocol')) | |
| gui_brainstorm('SetCurrentProtocol', iProtocol); | |
| end | |
| %% Prepare files | |
| % Find Recording files | |
| dataFiles = bst_process('CallProcess', 'process_select_files_data', [], [], ... | |
| 'subjectname', SubjectName, ... | |
| 'condition', ConditionName, ... | |
| 'tag', '', ... | |
| 'includebad', 0, ... | |
| 'includeintra', 0, ... | |
| 'includecommon', 0); | |
| % Find Surface Contrained source files | |
| constrFiles = bst_process('CallProcess', 'process_select_files_results', [], [], ... | |
| 'subjectname', SubjectName, ... | |
| 'condition', ConditionName, ... | |
| 'tag', '(surface)(constr)', ... | |
| 'includebad', 0, ... | |
| 'includeintra', 0, ... | |
| 'includecommon', 0); | |
| % Find Surface Uncontrained source files | |
| unconstrFiles = bst_process('CallProcess', 'process_select_files_results', [], [], ... | |
| 'subjectname', SubjectName, ... | |
| 'condition', ConditionName, ... | |
| 'tag', '(surface)(unconstr)', ... | |
| 'includebad', 0, ... | |
| 'includeintra', 0, ... | |
| 'includecommon', 0); | |
| %% Flattening | |
| % Apply PCA PCAi, PCAa and Norm flattening to Unconstrained and get files | |
| % Process: Unconstrained to flat map | |
| pcaFlatFiles = bst_process('CallProcess', 'process_source_flat', unconstrFiles, [], ... | |
| 'method', 'pca', ... % PCA: First mode of svd(x,y,z), maximizes retained power | |
| 'pcaedit', struct(... | |
| 'Method', 'pca', ... | |
| 'Baseline', [], ... | |
| 'DataTimeWindow', [], ... | |
| 'RemoveDcOffset', 'file')); | |
| % Process: Unconstrained to flat map | |
| pcaiFlatFiles = bst_process('CallProcess', 'process_source_flat', unconstrFiles, [], ... | |
| 'method', 'pca', ... % PCA: First mode of svd(x,y,z), maximizes retained power | |
| 'pcaedit', struct(... | |
| 'Method', 'pcai', ... | |
| 'Baseline', [], ... | |
| 'DataTimeWindow', [], ... | |
| 'RemoveDcOffset', 'file')); | |
| % Process: Unconstrained to flat map | |
| pcaaFlatFiles = bst_process('CallProcess', 'process_source_flat', unconstrFiles, [], ... | |
| 'method', 'pca', ... % PCA: First mode of svd(x,y,z), maximizes retained power | |
| 'pcaedit', struct(... | |
| 'Method', 'pcaa', ... | |
| 'Baseline', [], ... | |
| 'DataTimeWindow', [], ... | |
| 'RemoveDcOffset', 'file')); | |
| % Process: Unconstrained to flat map | |
| normFlatFiles = bst_process('CallProcess', 'process_source_flat', unconstrFiles, [], ... | |
| 'method', 'norm'); ... % Norm: sqrt(x^2+y^2+z^2) | |
| %% Sensor level | |
| % Process: Coherence 1xN [2021] | |
| sTmpFile = bst_process('CallProcess', 'process_cohere1_2021', dataFiles, [], ... | |
| 'timewindow', [], ... | |
| 'src_channel', emgChannel, ... | |
| 'dest_sensors', 'MEG', ... | |
| 'includebad', 0, ... | |
| 'removeevoked', 0, ... | |
| 'cohmeasure', 'mscohere', ... | |
| 'win_length', winLength, ... | |
| 'overlap', overlap, ... | |
| 'maxfreq', maxFreq, ... | |
| 'outputmode', 'avgcoh'); % Average cross-spectra of input files (one output file) | |
| % Process: Add tag | |
| bst_process('CallProcess', 'process_add_tag', sTmpFile, [], ... | |
| 'tag', '(sensors)', ... | |
| 'output', 1); % Add to file name | |
| %% Source level | |
| % For all types of source level files do: | |
| % 1. EMG vs Sources | |
| % 2. EMG vs Scouts (mean before) | |
| % 3. EMG vs Scouts (mean after ) | |
| % 4. EMG vs Scouts (pca before) | |
| % 5. EMG vs Scouts (pcai before) | |
| % 6. EMG vs Scouts (pcaa before) | |
| sourceTypesFiles = {'Constr', constrFiles; ... | |
| 'Unconstr', unconstrFiles; ... | |
| 'PCA flat', pcaFlatFiles; ... | |
| 'PCAi flat', pcaiFlatFiles; ... | |
| 'PCAa flat', pcaaFlatFiles; ... | |
| 'Norm flat', normFlatFiles}; | |
| for iType = 1 : size(sourceTypesFiles, 1) | |
| sourceType = ['(', sourceTypesFiles{iType, 1}, ')']; | |
| sourceFiles = sourceTypesFiles{iType, 2}; | |
| % 1. EMG vs Sources | |
| % Process: Coherence AxB [2021] | |
| sTmpFile = bst_process('CallProcess', 'process_cohere2_2021', dataFiles, sourceFiles, ... | |
| 'timewindow', [], ... | |
| 'src_channel', emgChannel, ... | |
| 'removeevoked', 0, ... | |
| 'cohmeasure', cohMeasure, ... | |
| 'win_length', winLength, ... | |
| 'overlap', overlap, ... | |
| 'maxfreq', maxFreq, ... | |
| 'outputmode', 'avgcoh'); % Average cross-spectra of input files (one output file) | |
| % Process: Add tag | |
| bst_process('CallProcess', 'process_add_tag', sTmpFile, [], ... | |
| 'tag', ['(sources)', sourceType], ... | |
| 'output', 1); % Add to file name | |
| % 2. EMG vs Scouts (mean before) | |
| % Process: Coherence AxB [2021] | |
| sTmpFile = bst_process('CallProcess', 'process_cohere2_2021', dataFiles, sourceFiles, ... | |
| 'timewindow', [], ... | |
| 'src_channel', emgChannel, ... | |
| 'dest_scouts', atlasScouts, ... | |
| 'flatten', 0, ... | |
| 'scouttime', 'before', ... % before | |
| 'scoutfunc', 'mean', ... % Mean | |
| 'scoutfuncaft', 'mean', ... % Mean | |
| 'pcaedit', struct(... | |
| 'Method', 'pca', ... | |
| 'Baseline', [], ... | |
| 'DataTimeWindow', [], ... | |
| 'RemoveDcOffset', 'file'), ... | |
| 'removeevoked', 0, ... | |
| 'cohmeasure', cohMeasure, ... | |
| 'win_length', winLength, ... | |
| 'overlap', overlap, ... | |
| 'maxfreq', maxFreq, ... | |
| 'outputmode', 'avgcoh'); % Average cross-spectra of input files (one output file) | |
| % Process: Add tag | |
| bst_process('CallProcess', 'process_add_tag', sTmpFile, [], ... | |
| 'tag', sourceType, ... | |
| 'output', 1); % Add to file name | |
| % 3. EMG vs Scouts (mean after ) | |
| sTmpFile = bst_process('CallProcess', 'process_cohere2_2021', dataFiles, sourceFiles, ... | |
| 'timewindow', [], ... | |
| 'src_channel', emgChannel, ... | |
| 'dest_scouts', atlasScouts, ... | |
| 'flatten', 0, ... | |
| 'scouttime', 'after', ... % after | |
| 'scoutfunc', 'mean', ... % Mean | |
| 'scoutfuncaft', 'mean', ... % Mean | |
| 'pcaedit', struct(... | |
| 'Method', 'pca', ... | |
| 'Baseline', [], ... | |
| 'DataTimeWindow', [], ... | |
| 'RemoveDcOffset', 'file'), ... | |
| 'removeevoked', 0, ... | |
| 'cohmeasure', cohMeasure, ... | |
| 'win_length', winLength, ... | |
| 'overlap', overlap, ... | |
| 'maxfreq', maxFreq, ... | |
| 'outputmode', 'avgcoh'); % Average cross-spectra of input files (one output file) | |
| % Process: Add tag | |
| bst_process('CallProcess', 'process_add_tag', sTmpFile, [], ... | |
| 'tag', sourceType, ... | |
| 'output', 1); % Add to file name | |
| % 4. EMG vs Scouts (pca before) | |
| sTmpFile = bst_process('CallProcess', 'process_cohere2_2021', dataFiles, sourceFiles, ... | |
| 'timewindow', [], ... | |
| 'src_channel', emgChannel, ... | |
| 'dest_scouts', atlasScouts, ... | |
| 'flatten', 0, ... | |
| 'scouttime', 'before', ... % before | |
| 'scoutfunc', 'pca', ... % PCA | |
| 'scoutfuncaft', 'mean', ... % Mean | |
| 'pcaedit', struct(... | |
| 'Method', 'pca', ... | |
| 'Baseline', [], ... | |
| 'DataTimeWindow', [], ... | |
| 'RemoveDcOffset', 'file'), ... | |
| 'removeevoked', 0, ... | |
| 'cohmeasure', cohMeasure, ... | |
| 'win_length', winLength, ... | |
| 'overlap', overlap, ... | |
| 'maxfreq', maxFreq, ... | |
| 'outputmode', 'avgcoh'); % Average cross-spectra of input files (one output file) | |
| % Process: Add tag | |
| bst_process('CallProcess', 'process_add_tag', sTmpFile, [], ... | |
| 'tag', ['(pca)', sourceType], ... | |
| 'output', 1); % Add to file name | |
| % 5. EMG vs Scouts (pcai before) | |
| sTmpFile = bst_process('CallProcess', 'process_cohere2_2021', dataFiles, sourceFiles, ... | |
| 'timewindow', [], ... | |
| 'src_channel', emgChannel, ... | |
| 'dest_scouts', atlasScouts, ... | |
| 'flatten', 0, ... | |
| 'scouttime', 'before', ... % before | |
| 'scoutfunc', 'pca', ... % PCA | |
| 'scoutfuncaft', 'mean', ... % Mean | |
| 'pcaedit', struct(... | |
| 'Method', 'pcai', ... | |
| 'Baseline', [], ... | |
| 'DataTimeWindow', [], ... | |
| 'RemoveDcOffset', 'file'), ... | |
| 'removeevoked', 0, ... | |
| 'cohmeasure', cohMeasure, ... | |
| 'win_length', winLength, ... | |
| 'overlap', overlap, ... | |
| 'maxfreq', maxFreq, ... | |
| 'outputmode', 'avgcoh'); % Average cross-spectra of input files (one output file) | |
| % Process: Add tag | |
| bst_process('CallProcess', 'process_add_tag', sTmpFile, [], ... | |
| 'tag', ['(pcai)', sourceType], ... | |
| 'output', 1); % Add to file name | |
| % 6. EMG vs Scouts (pcaa before) | |
| sTmpFile = bst_process('CallProcess', 'process_cohere2_2021', dataFiles, sourceFiles, ... | |
| 'timewindow', [], ... | |
| 'src_channel', emgChannel, ... | |
| 'dest_scouts', atlasScouts, ... | |
| 'flatten', 0, ... | |
| 'scouttime', 'before', ... % before | |
| 'scoutfunc', 'pca', ... % Mean | |
| 'scoutfuncaft', 'mean', ... % Mean | |
| 'pcaedit', struct(... | |
| 'Method', 'pcaa', ... | |
| 'Baseline', [], ... | |
| 'DataTimeWindow', [], ... | |
| 'RemoveDcOffset', 'file'), ... | |
| 'removeevoked', 0, ... | |
| 'cohmeasure', cohMeasure, ... | |
| 'win_length', winLength, ... | |
| 'overlap', overlap, ... | |
| 'maxfreq', maxFreq, ... | |
| 'outputmode', 'avgcoh'); % Average cross-spectra of input files (one output file) | |
| % Process: Add tag | |
| bst_process('CallProcess', 'process_add_tag', sTmpFile, [], ... | |
| 'tag', ['(pcaa)', sourceType], ... | |
| 'output', 1); % Add to file name | |
| end | |
| % Process: Send report by email | |
| bst_process('CallProcess', 'process_report_email', [], [], ... | |
| 'username', '', ... | |
| 'cc', '', ... | |
| 'subject', ['Process completed! Started at: ', startTime], ... | |
| 'reportfile', 'current', ... | |
| 'full', 1); | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment