make percentile threshold a parameter for tCompCor

Author: Shoshana Berleant

nipype/algorithms/compcor.py

@@ -102,8 +102,20 @@ def _add_extras(self, components, extra_regressors):
             regressors = np.genfromtxt(self.inputs.extra_regressors)
             return np.hstack((components, regressors))
 
+class TCompCorInputSpec(CompCorInputSpec):
+    # and all the fields in CompCorInputSpec
+    percentile_threshold = traits.Range(low=0., high=1., value=.02,
+                                        exclude_low=True, exclude_high=True,
+                                        usedefault=True, desc='the percentile '
+                                        'used to select highest-variance '
+                                        'voxels. By default the 2% of voxels '
+                                        'with the highest variance are used.')
+
 class TCompCor(CompCor):
 
+    input_spec = TCompCorInputSpec
+    output_spec = CompCorOutputSpec
+
     def _run_interface(self, runtime):
         imgseries = nb.load(self.inputs.realigned_file).get_data()
         time_voxels = imgseries.T
@@ -119,16 +131,16 @@ def _run_interface(self, runtime):
         tSTD = self._compute_tSTD(time_voxels, 0)
         sortSTD = np.sort(tSTD, axis=None) # flattened sorted matrix
 
-        # "... and retained a pre-specified upper fraction of the sorted voxels
-        # within each slice ... we chose a 2% threshold"
-        threshold = sortSTD[int(num_voxels * .98)]
-        mask = tSTD >= threshold
+        # use percentile_threshold to pick voxels
+        threshold_index = int(num_voxels * (1. - self.inputs.percentile_threshold))
+        threshold_std = sortSTD[threshold_index]
+        mask = tSTD >= threshold_std
         mask = mask.astype(int)
 
         # save mask
         mask_file = 'mask.nii'
         nb.nifti1.save(nb.Nifti1Image(mask, np.eye(4)), mask_file)
-        self.inputs.mask_file = 'mask.nii'
+        self.inputs.mask_file = mask_file
 
         super(TCompCor, self)._run_interface(runtime)
         return runtime

nipype/algorithms/tests/test_compcor.py

@@ -65,6 +65,14 @@ def test_tcompcor(self):
                                   ['-0.6048328569'], ['0.2134704201'],
                                   ['-0.0355784033']])
 
+    def test_tcompcor_with_percentile(self):
+        ccinterface = TCompCor(realigned_file=self.realigned_file, percentile_threshold=0.2)
+        ccinterface.run()
+
+        mask = nb.load('mask.nii').get_data()
+        num_nonmasked_voxels = np.count_nonzero(mask)
+        assert_equal(num_nonmasked_voxels, 2)
+
     def run_cc(self, ccinterface, expected_components):
         # run
         ccresult = ccinterface.run()