ENH: ANTs JointFusion with nosetests

CHANGES

@@ -22,6 +22,7 @@ Next release
 * FIX: OpenfMRI support and FSL 5.0.7 changes (https://github.com/nipy/nipype/pull/1006)
 * FIX: Output prefix in SPM Normalize with modulation (https://github.com/nipy/nipype/pull/1023)
 * ENH: Usability improvements in cluster environments (https://github.com/nipy/nipype/pull/1025)
+* ENH: ANTs JointFusion() (https://github.com/nipy/nipype/pull/1042)
 
 Release 0.10.0 (October 10, 2014)
 ============

nipype/interfaces/ants/__init__.py

@@ -11,7 +11,7 @@
 
 
 # Segmentation Programs
-from .segmentation import Atropos, LaplacianThickness, N4BiasFieldCorrection
+from .segmentation import Atropos, LaplacianThickness, N4BiasFieldCorrection, JointFusion
 
 # Visualization Programs
 from .visualization import ConvertScalarImageToRGB, CreateTiledMosaic

nipype/interfaces/ants/base.py

@@ -68,6 +68,11 @@ def _num_threads_update(self):
             self.inputs.environ.update({PREFERED_ITKv4_THREAD_LIMIT_VARIABLE:
                                             '%s' % self.inputs.num_threads})
 
+    @staticmethod
+    def _format_xarray(val):
+        """ Convenience method for converting input arrays [1,2,3] to commandline format '1x2x3' """
+        return 'x'.join([str(x) for x in val])
+
     @classmethod
     def set_default_num_threads(cls, num_threads):
         """Set the default number of threads for ITK calls

nipype/interfaces/ants/segmentation.py

@@ -570,3 +570,87 @@ def _list_outputs(self):
                                                                'subjectToTemplateLogJacobian.' +
                                                                self.inputs.image_suffix)
         return outputs
+
+class JointFusionInputSpec(ANTSCommandInputSpec):
+    dimension = traits.Enum(3, 2, 4, argstr='%d', position=0, usedefault=True, mandatory=True,
+                            desc='image dimension (2, 3, or 4)')
+    modalities = traits.Int(argstr='%d', position=1, mandatory=True, desc='Number of modalities or features')
+    warped_intensity_images = InputMultiPath(File(exists=True), argstr="-g %s...", mandatory=True, desc='Warped atlas images')
+    target_image = InputMultiPath(File(exists=True), argstr='-tg %s...', mandatory=True, desc='Target image(s)')
+    warped_label_images = InputMultiPath(File(exists=True), argstr="-l %s...", mandatory=True, desc='Warped atlas segmentations')
+    method = traits.Str(default='Joint', argstr='-m %s', usedefault=True, desc='Select voting method. Options: Joint (Joint Label Fusion). May be followed by optional parameters in brackets, e.g., -m Joint[0.1,2]')
+    alpha = traits.Float(default=0.1, usedefault=True, requires=['method'], desc='Regularization term added to matrix Mx for inverse')
+    beta = traits.Int(default=2, usedefault=True, requires=['method'], desc='Exponent for mapping intensity difference to joint error')
+    output_label_image = File(argstr='%s', mandatory=True, position=-1, desc='Output fusion label map image')
+    patch_radius = traits.ListInt(minlen=3, maxlen=3, argstr='-rp %s', desc='Patch radius for similarity measures, scalar or vector. Default: 2x2x2')
+    search_radius = traits.ListInt(minlen=3, maxlen=3, argstr='-rs %s', desc='Local search radius. Default: 3x3x3')
+    exclusion_region = File(exists=True, argstr='-x %s', desc='Specify an exclusion region for the given label.')
+    output_posteriors_name_template = traits.Str('POSTERIOR_%02d.nii.gz', argstr='-p %s',
+                                                 desc="Save the posterior maps (probability that each voxel belongs to each " +\
+                                                 "label) as images. The number of images saved equals the number of labels. " +\
+                                                 "The filename pattern must be in C printf format, e.g. posterior%04d.nii.gz")
+    output_voting_weights_name_template = traits.Str('WEIGHTED_%04d.nii.gz', argstr='-w %s', desc="Save the voting weights as " +\
+                                                     "images. The number of images saved equals the number of atlases. The " +\
+                                                     "filename pattern must be in C printf format, e.g. weight%04d.nii.gz")
+    atlas_group_id = traits.ListInt(argstr='-gp %d...', desc='Assign a group ID for each atlas')
+    atlas_group_weights = traits.ListInt(argstr='-gpw %d...', desc='Assign the voting weights to each atlas group')
+
+
+class JointFusionOutputSpec(TraitedSpec):
+    output_label_image = File(exists=True)
+    # TODO: optional outputs - output_posteriors, output_voting_weights
+
+
+class JointFusion(ANTSCommand):
+    """
+    Examples
+    --------
+
+    >>> from nipype.interfaces.ants import JointFusion
+    >>> at = JointFusion()
+    >>> at.inputs.dimension = 3
+    >>> at.inputs.modalities = 1
+    >>> at.inputs.method = 'Joint[0.1,2]'
+    >>> at.inputs.output_label_image ='fusion_labelimage_output.nii'
+    >>> at.inputs.warped_intensity_images = ['im1.nii',
+    ...                                      'im2.nii']
+    >>> at.inputs.warped_label_images = ['segmentation0.nii.gz',
+    ...                                  'segmentation1.nii.gz']
+    >>> at.inputs.target_image = 'T1.nii'
+    >>> at.inputs.patch_radius = [3,2,1]
+    >>> at.inputs.search_radius = [1,2,3]
+    >>> at.cmdline
+    'jointfusion 3 1 -m Joint[0.1,2] -rp 3x2x1 -rs 1x2x3 -tg T1.nii -g im1.nii -g im2.nii -l segmentation0.nii.gz -l segmentation1.nii.gz fusion_labelimage_output.nii'
+
+    Alternately, you can specify the voting method and parameters more 'Pythonically':
+
+    >>> at.inputs.method = 'Joint'
+    >>> at.inputs.alpha = 0.5
+    >>> at.inputs.beta = 1
+    >>> at.cmdline
+    'jointfusion 3 1 -m Joint[0.5,1] -rp 3x2x1 -rs 1x2x3 -tg T1.nii -g im1.nii -g im2.nii -l segmentation0.nii.gz -l segmentation1.nii.gz fusion_labelimage_output.nii'
+    """
+    input_spec = JointFusionInputSpec
+    output_spec = JointFusionOutputSpec
+    _cmd = 'jointfusion'
+
+    def _format_arg(self, opt, spec, val):
+        if opt == 'method':
+            if '[' in val:
+                retval = '-m {0}'.format(val)
+            else:
+                retval = '-m {0}[{1},{2}]'.format(self.inputs.method, self.inputs.alpha, self.inputs.beta)
+        elif opt == 'patch_radius':
+            retval = '-rp {0}'.format(self._format_xarray(val))
+        elif opt == 'search_radius':
+            retval = '-rs {0}'.format(self._format_xarray(val))
+        else:
+            if opt == 'warped_intensity_images':
+                assert len(val) == len(self.inputs.warped_label_images), "Number of intensity images and label maps must be the same"
+            return super(ANTSCommand, self)._format_arg(opt, spec, val)
+        return retval
+
+    def _list_outputs(self):
+        outputs = self._outputs().get()
+        outputs['output_label_image'] = os.path.abspath(self.inputs.output_label_image)
+        return outputs

nipype/interfaces/ants/tests/test_JointFusion.py

@@ -0,0 +1,77 @@
+from nose import with_setup
+from nipype.testing import assert_equal, assert_raises
+from nipype.interfaces.base import InputMultiPath
+from traits.trait_errors import TraitError
+from nipype.interfaces.ants import JointFusion
+
+
+def test_JointFusion_dimension():
+    at = JointFusion()
+    set_dimension = lambda d: setattr(at.inputs, 'dimension', int(d))
+    for d in range(2, 5):
+        set_dimension(d)
+        yield assert_equal, at.inputs.dimension, int(d)
+    for d in [0, 1, 6, 7]:
+        yield assert_raises, TraitError, set_dimension, d
+
+
+def test_JointFusion_modalities():
+    at = JointFusion()
+    set_modalities = lambda m: setattr(at.inputs, 'modalities', int(m))
+    for m in range(1, 5):
+        set_modalities(m)
+        yield assert_equal, at.inputs.modalities, int(m)
+
+
+def test_JointFusion_method():
+    at = JointFusion()
+    set_method = lambda a, b: setattr(at.inputs, 'method', 'Joint[%.1f,%d]'.format(a, b))
+    for a in range(10):
+        _a = a / 10.0
+        for b in range(10):
+            set_method(_a, b)
+            # set directly
+            yield assert_equal, at.inputs.method, 'Joint[%.1f,%d]'.format(_a, b)
+            aprime = _a + 0.1
+            bprime = b + 1
+            at.inputs.alpha = aprime
+            at.inputs.beta = bprime
+            # set with alpha/beta
+            yield assert_equal, at.inputs.method, 'Joint[%.1f,%d]'.format(aprime, bprime)
+
+
+def test_JointFusion_radius():
+    at = JointFusion()
+    set_radius = lambda attr,x,y,z: setattr(at.inputs, attr, [x, y, z])
+    for attr in ['patch_radius', 'search_radius']:
+        for x in range(5):
+            set_radius(attr, x, x + 1, x**x)
+            yield assert_equal, at._format_arg(attr, None, getattr(at.inputs, attr))[4:], '{0}x{1}x{2}'.format(x, x + 1, x**x)
+
+
+def test_JointFusion_cmd():
+    at = JointFusion()
+    at.inputs.dimension = 3
+    at.inputs.modalities = 1
+    at.inputs.method = 'Joint[0.1,2]'
+    at.inputs.output_label_image ='fusion_labelimage_output.nii'
+    at.inputs.warped_intensity_images = ['im1.nii',
+                                         'im2.nii']
+    at.inputs.warped_label_images = ['segmentation0.nii.gz',
+                                     'segmentation1.nii.gz']
+    at.inputs.target_image = 'T1.nii'
+    at.inputs.patch_radius = [3,2,1]
+    at.inputs.search_radius = [1,2,3]
+    yield assert_equal, at.cmdline, 'jointfusion 3 1' + \
+      ' -m Joint[0.1,2]' + \
+      ' -rp 3x2x1' + \
+      ' -rs 1x2x3' + \
+      ' -tg T1.nii' + \
+      ' -g im1.nii' + \
+      ' -g im2.nii' + \
+      ' -l segmentation0.nii.gz' + \
+      ' -l segmentation1.nii.gz' + \
+      ' fusion_labelimage_output.nii',
+    "Command line generation failure", True
+    # setting intensity or labels with unequal lengths raises error
+    yield assert_raises, AssertionError, at._format_arg, 'warped_intensity_images', InputMultiPath, ['im1.nii', 'im2.nii', 'im3.nii']