add efficiency demo

cpp-lln-lab · Feb 16, 2024 · f84be7e · f84be7e
1 parent 39a55c3
commit f84be7e
Show file tree

Hide file tree

Showing 5 changed files with 337 additions and 269 deletions.
diff --git a/demos/efficiency/demo_efficiency.m b/demos/efficiency/demo_efficiency.m
@@ -1,34 +1,61 @@
 % (C) Copyright 2023 bidspm developers
 
+% Example to compute efficiency from
+% https://www.nature.com/articles/s41598-024-52967-8
+%
+% Related to this thread:
+% https://twitter.com/m_wall/status/1757311725777387739
+
 clear;
 clc;
+close all;
 
 addpath(fullfile(pwd, '..', '..'));
 bidspm();
 
+conditions = {'alcohol', 'neutral', 'negative', 'positive'};
+for i = 1:numel(conditions)
+  design_matrix{i} = ['trial_type.' conditions{i}];
+end
+
 %% create stats model JSON
 bm = BidsModel();
 NodeIdx = 1;
-bm.Nodes{NodeIdx}.Model.X = {'trial_type.cdt_A', 'trial_type.cdt_B'};
-bm.Nodes{NodeIdx}.DummyContrasts = struct('type', 't', ...
-                                          'Contrasts', {{'trial_type.cdt_A', 'trial_type.cdt_B'}});
 
-contrast = struct('type', 't', ...
-                  'Name', 'A_gt_B', ...
-                  'Weights', [1, -1], ...
-                  'ConditionList', {{'trial_type.cdt_A', 'trial_type.cdt_B'}});
+% high pass filter
+bm.Nodes{NodeIdx}.Model.Options.HighPassFilterCutoffHz = 1 / 128;
 
-bm.Nodes{NodeIdx}.Contrasts{1} = contrast;
+for i = 1:numel(conditions)
+  bm.Nodes{NodeIdx}.Model.X = design_matrix;
+end
+
+% contrast for each condition against baseline
+bm.Nodes{NodeIdx}.DummyContrasts = struct('type', 't', ...
+                                          'Contrasts', {design_matrix});
+
+% contrast of interests
+% alcohol > neural
+% negative > neutral
+% positive > neutral
+conditions_to_check = {'alcohol', 'negative', 'positive'};
+for i = 1:numel(conditions_to_check)
+  contrast = struct('Test', 't', ...
+                    'Name', [conditions_to_check{i} '_gt_neutral'], ...
+                    'Weights', [1, -1], ...
+                    'ConditionList', {{['trial_type.' conditions_to_check{i}], 'trial_type.neutral'}});
+  bm.Nodes{NodeIdx}.Contrasts{i} = contrast;
+end
 
 bm.write('smdl.json');
 
 %% create events TSV file
-conditions = {'cdt_A', 'cdt_B'};
-IBI = 5;
-ISI = 0.1;
-stimDuration = 1.5;
-stimPerBlock = 12;
-nbBlocks = 10;
+IBI = 25;
+ISI = 0;
+stimDuration = 4;
+stimPerBlock = 40;
+nbBlocks = 1;
+
+TR = 1;
 
 trial_type = {};
 onset = [];
@@ -52,8 +79,15 @@
 
 bids.util.tsvwrite('events.tsv', tsv);
 
-opt.TR = 2;
+opt.TR = TR;
 
 opt.model.file = fullfile(pwd, 'smdl.json');
 
-e = computeDesignEfficiency(fullfile(pwd, 'events.tsv'), opt);
+[e, X] = computeDesignEfficiency(fullfile(pwd, 'events.tsv'), opt);
+
+%%
+% saveas(gcf(), 'design_matrix.png');
+% close(gcf());
+%
+% saveas(gcf(), 'events.png');
+% close(gcf());