Fixes for HBN data to be compatible

fabriciocravo · fabriciocravo · commit d44e50e8fcaa · 2026-03-02T16:24:21.000-03:00
diff --git a/.DS_Store b/.DS_Store
diff --git a/parfor_stat_tools/pf_repetition_loop.m b/parfor_stat_tools/pf_repetition_loop.m
@@ -80,7 +80,7 @@
             end
         
             % Skip method if no submethods need computation
-            if ~any(struct2array(submethods_struct))
+            if ~any(cell2mat(struct2cell(submethods_struct)))
                 continue;
             end
         
diff --git a/power_calculator_tools/atlas_data_set_map.m b/power_calculator_tools/atlas_data_set_map.m
@@ -58,7 +58,9 @@ case startsWith(data_set_choice, 'test') && endsWith(data_set_choice, 'act')
                 atlas_file = NaN;
 
             otherwise
-                error('No atlas file found for dataset: %s', data_set_choice);
+                % Default to shen atlas now
+                atlas_file = './atlas_storage/map268_subnetwork.mat';
+                % error('No atlas file found for dataset: %s', data_set_choice);
         end
     end
 end
diff --git a/power_calculator_tools/cap_t2_test_subs_when_data_limited.m b/power_calculator_tools/cap_t2_test_subs_when_data_limited.m
@@ -0,0 +1,16 @@
+function subset_cell = cap_t2_test_subs_when_data_limited(RP)
+
+    subset_cell = {};
+    for i = 1:numel(RP.list_of_nsubset)
+        n = RP.list_of_nsubset{i};
+
+        if n >= RP.n_subs_1 || n >= (RP.n_subs - RP.n_subs_1)
+            fprintf(['Subset size %d exceeds minority group (%d vs %d). ' ...
+                'Removing from analysis.\n'], ...
+                n, RP.n_subs_1, RP.n_subs - RP.n_subs_1);
+        else
+            subset_cell{end+1} = n;
+        end
+    end
+
+end
diff --git a/power_calculator_tools/draw_repetition_ids.m b/power_calculator_tools/draw_repetition_ids.m
@@ -57,7 +57,7 @@
                 ids = [ids; ids + RP.n_subs]';
     
             case 't2'
-                
+
                 ids_1 = randperm(RP.n_subs_1, RP.n_subs_subset);
                 ids_2 = randperm(RP.n_subs - RP.n_subs_1, RP.n_subs_subset) + RP.n_subs_1;
     
diff --git a/power_calculator_tools/extract_atlas_related_parameters.m b/power_calculator_tools/extract_atlas_related_parameters.m
@@ -43,7 +43,7 @@
     trilmask_net = NaN;
     edge_groups = [];
     n_networks = 0;
-
+       
     % Only apply atlas to the network-based stats
     if ~isnan(RP.atlas_file)
 
diff --git a/power_calculator_tools/get_index_matching_score.m b/power_calculator_tools/get_index_matching_score.m
@@ -0,0 +1,11 @@
+function [index_cond_1, index_cond_2] = get_index_matching_score(TestData_score, test_score_set)
+    
+    % Type cast to str to handle weird datatypes
+    score_str = string(TestData_score);
+    set_str = string(test_score_set);
+        
+    % Find indexes matching the condition
+    index_cond_1 = score_str == set_str(1);
+    index_cond_2 = score_str == set_str(2);
+
+end
diff --git a/power_calculator_tools/infer_test_from_data.m b/power_calculator_tools/infer_test_from_data.m
@@ -55,19 +55,17 @@
 
         test_type_origin = 'score_cond';
         % if two unique entries in score -> t (paired) or t2
-
-        % Are the subids always in outcomes? For the hpc one, we have
-        % NaN and it is likely a t2 test
-        index_cond_1 = strcmp(TestData.score, test_score_set{1});
-        index_cond_2 = strcmp(TestData.score, test_score_set{2});
+        
+        [index_cond_1, index_cond_2] = get_index_matching_score(TestData.score, ...
+            test_score_set);
 
         sub_ids_cond1 = BrainData.(TestData.reference_condition).sub_ids(index_cond_1);
         sub_ids_cond2 = BrainData.(TestData.reference_condition).sub_ids(index_cond_2);
         
         %% TODO: Divided by the two - focus on group sizes
         n_equal = numel(intersect(sort(sub_ids_cond1), sort(sub_ids_cond2)));
         n_unique = numel(setxor(sub_ids_cond1, sub_ids_cond2));
- 
+
         if n_equal >= n_unique
             test_type = 't';
         else
diff --git a/power_calculator_tools/rep_cal_function.m b/power_calculator_tools/rep_cal_function.m
@@ -106,7 +106,7 @@ function rep_cal_function(Params)
         % However, only the r test type will use this test
         switch test_type_origin
             
-            % Here: RP.test_name, RP.n_subs_1, RP.n_subs_2, RP.n_subs
+            % Here: RP.test_name, RP.n_subs_1, RP.n_subs_2, RP.n_subs 
             case 'score_cond'
                 [X, Y , RP] = subs_data_from_score_condition(RP, OutcomeData.(t), BrainData, t);
            
@@ -117,6 +117,10 @@ function rep_cal_function(Params)
                 error('Test type origin not found')
 
         end  
+        
+        if strcmp(RP.test_type, 't2')
+            RP.list_of_nsubset = cap_t2_test_subs_when_data_limited(RP);
+        end
 
         % Sets ground truth parameters for gt calculation
         % Future note, avoid this, for clarity, make sure each function returns exactly
diff --git a/power_calculator_tools/subs_data_from_score_condition.m b/power_calculator_tools/subs_data_from_score_condition.m
@@ -43,10 +43,10 @@
         case 't2'
            
             ref_cond = TestData.reference_condition;
-    
-            index_cond_1 = strcmp(TestData.score, test_score_set{1});
-            index_cond_2 = strcmp(TestData.score, test_score_set{2});
             
+            [index_cond_1, index_cond_2] = get_index_matching_score(TestData.score, ...
+                test_score_set);
+
             % I call c1 one rest and c2 task for consitency - even though
             % they might not be task and rest
             RP.sub_ids_rest = TestData.sub_ids(index_cond_1);
diff --git a/setparams.m b/setparams.m
@@ -29,11 +29,11 @@
 
 % Datasets - Commented for easy use
 % Params.data_dir = './data/s_abcd_fc_rosenblatt.mat';
-Params.data_dir = './data/s_hcp_fc_noble_tasks.mat';
-% Params.data_dir = './data/s_hcp_act_noble_1.mat';
+% Params.data_dir = './data/s_hcp_fc_noble_tasks.mat';
+Params.data_dir = './data/s_hbn_fc_ye.mat';
 
 % Name of the output directory 
-Params.output = 're_tfce_speed_comp';
+% Params.output = nan;
 
 % Save specifications - if NaN output becomes dataset file name
 Params.save_directory = './power_calculator_results/';
@@ -60,8 +60,8 @@
 %%% Resampling parameters %%%
 Params.parallel = true; % run stuff sequentially or in parallel
 Params.n_workers = 5; % num parallel workers for parfor, best if # workers = # cores
-Params.n_repetitions = 5;  % 500 recommended
-Params.batch_size = 5;
+Params.n_repetitions = 100;  % 500 recommended
+Params.batch_size = 10;
  
 
 %% Skip some tests - change ranges or the function
@@ -70,20 +70,20 @@
 
 
 %% List of subjects per subset
-Params.list_of_nsubset = {20, 80, 200}; % To change this, add more when necessary
+Params.list_of_nsubset = {20, 40, 80, 120, 200}; % To change this, add more when necessary
                     % size of subset is full group size (N=n*2 for two sample t-test or N=n for one-sample)
 
                             % Current model (see above design matrix) only designed for t-test
 Params.force_permute = true;               
-Params.n_perms = 100;               % recommend n_perms=5000 to appreciably reduce uncertainty of p-value estimation (https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/Randomise/Theory)
+Params.n_perms = 1000;               % recommend n_perms=5000 to appreciably reduce uncertainty of p-value estimation (https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/Randomise/Theory)
 Params.tthresh_first_level = 3.1;    % t=3.1 corresponds with p=0.005-0.001 (DOF=10-1000)
                             % Only used if cluster_stat_type='Size'
 Params.pthresh_second_level = 0.05;  % FWER or FDR rate 
 Params.tpr_dthresh = 0; % Threshold for true positives vs negatives
 Params.save_significance_thresh = 0.15;
-%Params.all_cluster_stat_types = {'Parametric', 'Size_cpp', 'Fast_TFCE_cpp', 'Constrained_cpp', 'Omnibus_cNBS'};
-Params.all_cluster_stat_types = {'IC_TFCE_FC_cpp_dh1', 'IC_TFCE_FC_cpp_dh5', 'IC_TFCE_FC_cpp_dh10',...
-    'IC_TFCE_FC_cpp_dh25', 'TFCE_cpp_dh1', 'TFCE_cpp_dh5','TFCE_cpp_dh10', 'TFCE_cpp_dh25', 'Exact_FC_TFCE_cpp'};
+Params.all_cluster_stat_types = {'Parametric', 'Size_cpp', 'Fast_TFCE_cpp', 'Constrained_cpp', 'Omnibus_cNBS'};
+%Params.all_cluster_stat_types = {'IC_TFCE_FC_cpp_dh1', 'IC_TFCE_FC_cpp_dh5', 'IC_TFCE_FC_cpp_dh10',...
+%    'IC_TFCE_FC_cpp_dh25', 'TFCE_cpp_dh1', 'TFCE_cpp_dh5','TFCE_cpp_dh10', 'TFCE_cpp_dh25', 'Exact_FC_TFCE_cpp'};
 
 Params.all_submethods = {'FWER', 'FDR'};
 
