@@ -586,6 +586,16 @@ class ThresholdInputSpec(SPMCommandInputSpec):
586586            "set to p-value)" 
587587        ),
588588    )
589+     use_vox_fdr_correction  =  traits .Bool (
590+         False ,
591+         usedefault = True ,
592+         desc = (
593+             "whether to use voxel-based FDR " 
594+             "correction for initial threshold " 
595+             "(height_threshold_type has to be " 
596+             "set to q-value)" 
597+         ),
598+     )
589599    use_topo_fdr  =  traits .Bool (
590600        True ,
591601        usedefault = True ,
@@ -661,8 +671,16 @@ def _gen_pre_topo_map_filename(self):
661671    def  _make_matlab_command (self , _ ):
662672        script  =  "con_index = %d;\n "  %  self .inputs .contrast_index 
663673        script  +=  "cluster_forming_thr = %f;\n "  %  self .inputs .height_threshold 
664-         if  self .inputs .use_fwe_correction :
674+ 
675+         if  self .inputs .use_fwe_correction  and  self .inputs .use_vox_fdr_correction :
676+             raise  ValueError (
677+                 "'use_fwe_correction' and 'use_vox_fdr_correction' can't both be True" 
678+             )
679+ 
680+         if  self .inputs .use_fwe_correction  and  not  self .inputs .use_vox_fdr_correction :
665681            script  +=  "thresDesc  = 'FWE';\n " 
682+         elif  self .inputs .use_vox_fdr_correction  and  not  self .inputs .use_fwe_correction :
683+             script  +=  "thresDesc  = 'FDR';\n " 
666684        else :
667685            script  +=  "thresDesc  = 'none';\n " 
668686
@@ -687,6 +705,8 @@ def _make_matlab_command(self, _):
687705FWHM  = SPM.xVol.FWHM; 
688706df = [SPM.xCon(con_index).eidf SPM.xX.erdf]; 
689707STAT = SPM.xCon(con_index).STAT; 
708+ VspmSv   = cat(1,SPM.xCon(con_index).Vspm); 
709+ 
690710R = SPM.xVol.R; 
691711S = SPM.xVol.S; 
692712n = 1; 
@@ -695,6 +715,9 @@ def _make_matlab_command(self, _):
695715    case 'FWE' 
696716        cluster_forming_thr = spm_uc(cluster_forming_thr,df,STAT,R,n,S); 
697717
718+     case 'FDR' 
719+         cluster_forming_thr = spm_uc_FDR(cluster_forming_thr,df,STAT,n,VspmSv,0); 
720+ 
698721    case 'none' 
699722        if strcmp(height_threshold_type, 'p-value') 
700723            cluster_forming_thr = spm_u(cluster_forming_thr^(1/n),df,STAT); 
0 commit comments