Merge pull request #1 from neuroprismlab/state_of_power

Preparing for discovery transfer

Author: fabriciocravo

NBS_benchmarking/setup_benchmarking.m

@@ -21,6 +21,8 @@
     RepParams.n_workers = RepParams.test_n_workers;
 end
 
+% Get the level - edge-level, network-level, or whole brain level
+RepParams.stat_level = set_statistic_level(RepParams.cluster_stat_type);
 
 %% Adjust NBS Method According to Stat Type
 % special setup for nonparametric FDR (part of classic NBS toolbox) and parametric FDR and Bonferroni (newly added here)  since can't run it using cluster_stat_type

atlas_storage/test_hcp_fc_atlas.mat

@@ -62,8 +62,7 @@
 Params.n_subs_subset = 40;   % 40 | 80 | 120
 
 %% List of subjects per subset
-Params.list_of_nsubset = {20, 40, 80, 120}; % To change this, add more when necessary
-Params.list_of_nsubset = {40}; % Redefined because I need a base with all 40 first
+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)
 
 
@@ -111,7 +110,7 @@
 %%%%% DEVELOPERS ONLY %%%%%
 % Use a small subset of permutations for faster development -- inappropriate for inference
 
-Params.testing = true;
+Params.testing = false;
 Params.test_n_perms = '100';
 Params.test_n_repetitions = 1;
 Params.test_n_workers = 10;

config_files/setparams.m

@@ -12,8 +12,8 @@ function rep_cal_function(Params)
 
     %% Parallel Workers 
     % Uncoment the disp line if setup is commented out - as reminder 
-    % setup_parallel_workers(Params.parallel, Params.n_workers);
-    disp('Debugging: Setup parallel workers deactived')
+    setup_parallel_workers(Params.parallel, Params.n_workers);
+    % disp('Debugging: Setup parallel workers deactived')
 
     OutcomeData = Dataset.outcome;
     BrainData = Dataset.brain_data;

power_calculator_tools/rep_cal_function.m

@@ -122,11 +122,11 @@ function run_benchmarking(RP, Y)
                 %% Run NBS repetitions
 
                 % Be careful with this parfor commented lol
-                % parfor (i_rep=1: RP.n_repetitions)
+                parfor (i_rep=1: RP.n_repetitions)
 
                 % Do not erase the comment bellow, it is used for debugging
                 % when running not in parallel 
-disp('Parof');  for i_rep = 1:RP.n_repetitions
+% disp('Parof');  for i_rep = 1:RP.n_repetitions
 
                     % Encapsulation of the most computationally intensive loop
                     [FWER_rep, edge_stats_all_rep, pvals_all_rep, cluster_stats_all_rep, ...
@@ -141,13 +141,21 @@ function run_benchmarking(RP, Y)
                     edge_stats_all_neg(:,i_rep) = edge_stats_all_neg_rep;
                     pvals_all_neg(:,i_rep) = pvals_all_neg_rep;
 
-                    if ~isempty(cluster_stats_all_rep) || ~isempty(cluster_stats_all_neg_rep)
+                    % If it is an edge method, no need for the
+                    % network stats  
+                    if ~strcmp(RP.stat_level, 'edge')
                         cluster_stats_all(:,:,i_rep) = cluster_stats_all_rep;
                         cluster_stats_all_neg(:,:,i_rep) = cluster_stats_all_neg_rep;
                     end
 
                 end
 
+                % An NaN for network-level in the edge case
+                if strcmp(RP.stat_level, 'edge')
+                    cluster_stats_all = NaN;
+                    cluster_stats_all_neg = NaN;
+                end
+                
                 if contains(UI.statistic_type.ui,'Constrained') || strcmp(UI.statistic_type.ui,'SEA') ...
                     || strcmp(UI.statistic_type.ui,'Omnibus')
                     cluster_stats_all = squeeze(cluster_stats_all);