Merge pull request #2167 from mgxd/enh/calcmedian

enh: calcmedian interface

Author: satra

examples/fmri_ants_openfmri.py

@@ -26,7 +26,7 @@
 import nipype.pipeline.engine as pe
 import nipype.algorithms.modelgen as model
 import nipype.algorithms.rapidart as ra
-from nipype.algorithms.misc import TSNR
+from nipype.algorithms.misc import TSNR, CalculateMedian
 from nipype.interfaces.c3 import C3dAffineTool
 from nipype.interfaces import fsl, Function, ants, freesurfer as fs
 import nipype.interfaces.io as nio
@@ -55,33 +55,6 @@
     'from scipy.special import legendre'
 ]
 
-def median(in_files):
-    """Computes an average of the median of each realigned timeseries
-
-    Parameters
-    ----------
-
-    in_files: one or more realigned Nifti 4D time series
-
-    Returns
-    -------
-
-    out_file: a 3D Nifti file
-    """
-    average = None
-    for idx, filename in enumerate(filename_to_list(in_files)):
-        img = nb.load(filename, mmap=NUMPY_MMAP)
-        data = np.median(img.get_data(), axis=3)
-        if average is None:
-            average = data
-        else:
-            average = average + data
-    median_img = nb.Nifti1Image(average / float(idx + 1), img.affine,
-                                img.header)
-    filename = os.path.join(os.getcwd(), 'median.nii.gz')
-    median_img.to_filename(filename)
-    return filename
-
 
 def create_reg_workflow(name='registration'):
     """Create a FEAT preprocessing workflow together with freesurfer
@@ -818,11 +791,7 @@ def check_behav_list(behav, run_id, conds):
     wf.connect(preproc, "outputspec.realigned_files", tsnr, "in_file")
 
     # Compute the median image across runs
-    calc_median = Node(Function(input_names=['in_files'],
-                                output_names=['median_file'],
-                                function=median,
-                                imports=imports),
-                       name='median')
+    calc_median = Node(CalculateMedian(), name='median')
     wf.connect(tsnr, 'detrended_file', calc_median, 'in_files')
 
     """

examples/rsfmri_vol_surface_preprocessing.py

@@ -67,7 +67,7 @@
 # mlab.MatlabCommand.set_default_paths('/software/matlab/spm12')
 
 from nipype.algorithms.rapidart import ArtifactDetect
-from nipype.algorithms.misc import TSNR
+from nipype.algorithms.misc import TSNR, CalculateMedian
 from nipype.interfaces.utility import Rename, Merge, IdentityInterface
 from nipype.utils.filemanip import filename_to_list
 from nipype.interfaces.io import DataSink, FreeSurferSource
@@ -623,11 +623,7 @@ def create_workflow(files,
     wf.connect(slice_timing, 'timecorrected_files', tsnr, 'in_file')
 
     # Compute the median image across runs
-    calc_median = Node(Function(input_names=['in_files'],
-                                output_names=['median_file'],
-                                function=median,
-                                imports=imports),
-                       name='median')
+    calc_median = Node(CalculateMedian(), name='median')
     wf.connect(tsnr, 'detrended_file', calc_median, 'in_files')
 
     """Segment and Register

examples/rsfmri_vol_surface_preprocessing_nipy.py

@@ -65,7 +65,7 @@
 from nipype import Workflow, Node, MapNode
 
 from nipype.algorithms.rapidart import ArtifactDetect
-from nipype.algorithms.misc import TSNR
+from nipype.algorithms.misc import TSNR, CalculateMedian
 from nipype.algorithms.confounds import ACompCor
 from nipype.interfaces.utility import Rename, Merge, IdentityInterface
 from nipype.utils.filemanip import filename_to_list
@@ -556,11 +556,7 @@ def create_workflow(files,
     wf.connect(realign, "out_file", tsnr, "in_file")
 
     # Compute the median image across runs
-    calc_median = Node(Function(input_names=['in_files'],
-                                output_names=['median_file'],
-                                function=median,
-                                imports=imports),
-                       name='median')
+    calc_median = Node(CalculateMedian(), name='median')
     wf.connect(tsnr, 'detrended_file', calc_median, 'in_files')
 
     """Segment and Register

nipype/algorithms/misc.py

@@ -33,7 +33,7 @@
                                InputMultiPath, OutputMultiPath,
                                BaseInterfaceInputSpec, isdefined,
                                DynamicTraitedSpec, Undefined)
-from ..utils.filemanip import fname_presuffix, split_filename
+from ..utils.filemanip import fname_presuffix, split_filename, filename_to_list
 from ..utils import NUMPY_MMAP
 
 from . import confounds
@@ -1380,6 +1380,94 @@ def merge_rois(in_files, in_idxs, in_ref,
     return out_file
 
 
+class CalculateMedianInputSpec(BaseInterfaceInputSpec):
+    in_files = InputMultiPath(File(exists=True, mandatory=True,
+      desc="One or more realigned Nifti 4D timeseries"))
+    median_file = traits.Str(desc="Filename prefix to store median images")
+    median_per_file = traits.Bool(False, usedefault=True,
+      desc="Calculate a median file for each Nifti")
+
+class CalculateMedianOutputSpec(TraitedSpec):
+    median_files = OutputMultiPath(File(exists=True),
+      desc="One or more median images")
+
+class CalculateMedian(BaseInterface):
+    """
+    Computes an average of the median across one or more 4D Nifti timeseries
+
+    Example
+    -------
+
+    >>> from nipype.algorithms.misc import CalculateMedian
+    >>> mean = CalculateMedian()
+    >>> mean.inputs.in_files = 'functional.nii'
+    >>> mean.run() # doctest: +SKIP
+
+    """
+    input_spec = CalculateMedianInputSpec
+    output_spec = CalculateMedianOutputSpec
+
+    def __init__(self, *args, **kwargs):
+        super(CalculateMedian, self).__init__(*args, **kwargs)
+        self._median_files = []
+
+    def _gen_fname(self, suffix, idx=None, ext=None):
+        if idx:
+            in_file = self.inputs.in_files[idx]
+        else:
+            if isinstance(self.inputs.in_files, list):
+                in_file = self.inputs.in_files[0]
+            else:
+                in_file = self.inputs.in_files
+        fname, in_ext = op.splitext(op.basename(in_file))
+        if in_ext == '.gz':
+            fname, in_ext2 = op.splitext(fname)
+            in_ext = in_ext2 + in_ext
+        if ext is None:
+            ext = in_ext
+        if ext.startswith('.'):
+            ext = ext[1:]
+        if self.inputs.median_file:
+            outname = self.inputs.median_file
+        else:
+            outname = '{}_{}'.format(fname, suffix)
+        if idx:
+            outname += str(idx)
+        return op.abspath('{}.{}'.format(outname, ext))
+
+    def _run_interface(self, runtime):
+        total = None
+        self._median_files = []
+        for idx, fname in enumerate(filename_to_list(self.inputs.in_files)):
+            img = nb.load(fname, mmap=NUMPY_MMAP)
+            data = np.median(img.get_data(), axis=3)
+            if self.inputs.median_per_file:
+                self._median_files.append(self._write_nifti(img, data, idx))
+            else:
+                if total is None:
+                    total = data
+                else:
+                    total += data
+        if not self.inputs.median_per_file:
+            self._median_files.append(self._write_nifti(img, total, idx))
+        return runtime
+
+    def _list_outputs(self):
+        outputs = self._outputs().get()
+        outputs['median_files'] = self._median_files
+        return outputs
+
+    def _write_nifti(self, img, data, idx, suffix='median'):
+        if self.inputs.median_per_file:
+            median_img = nb.Nifti1Image(data, img.affine, img.header)
+            filename = self._gen_fname(suffix, idx=idx)
+        else:
+            median_img = nb.Nifti1Image(data/(idx+1), img.affine, img.header)
+            filename = self._gen_fname(suffix)
+        median_img.to_filename(filename)
+        return filename
+
+
 # Deprecated interfaces ------------------------------------------------------
 
 class Distance(nam.Distance):

nipype/algorithms/tests/test_auto_CalculateMedian.py

@@ -0,0 +1,29 @@
+# AUTO-GENERATED by tools/checkspecs.py - DO NOT EDIT
+from __future__ import unicode_literals
+from ..misc import CalculateMedian
+
+
+def test_CalculateMedian_inputs():
+    input_map = dict(ignore_exception=dict(nohash=True,
+    usedefault=True,
+    ),
+    in_files=dict(),
+    median_file=dict(),
+    median_per_file=dict(usedefault=True,
+    ),
+    )
+    inputs = CalculateMedian.input_spec()
+
+    for key, metadata in list(input_map.items()):
+        for metakey, value in list(metadata.items()):
+            assert getattr(inputs.traits()[key], metakey) == value
+
+
+def test_CalculateMedian_outputs():
+    output_map = dict(median_files=dict(),
+    )
+    outputs = CalculateMedian.output_spec()
+
+    for key, metadata in list(output_map.items()):
+        for metakey, value in list(metadata.items()):
+            assert getattr(outputs.traits()[key], metakey) == value

nipype/algorithms/tests/test_misc.py

@@ -9,6 +9,8 @@
 from nipype.algorithms import misc
 from nipype.utils.filemanip import fname_presuffix
 from nipype.testing.fixtures import create_analyze_pair_file_in_directory
+from nipype.utils import NUMPY_MMAP
+from nipype.testing import example_data
 
 
 def test_CreateNifti(create_analyze_pair_file_in_directory):
@@ -31,4 +33,16 @@ def test_CreateNifti(create_analyze_pair_file_in_directory):
     result = create_nifti.run()
 
     assert os.path.exists(result.outputs.nifti_file)
-    assert nb.load(result.outputs.nifti_file)
+    assert nb.load(result.outputs.nifti_file, mmap=NUMPY_MMAP)
+
+def test_CalculateMedian(create_analyze_pair_file_in_directory):
+
+    mean = misc.CalculateMedian()
+
+    with pytest.raises(TypeError): mean.run()
+
+    mean.inputs.in_files = example_data('ds003_sub-01_mc.nii.gz')
+    eg = mean.run()
+
+    assert os.path.exists(eg.outputs.median_files)
+    assert nb.load(eg.outputs.median_files, mmap=NUMPY_MMAP)