[ENH] Add a FramewiseDisplacement interface

CHANGES

@@ -1,6 +1,7 @@
 Upcoming release 0.13
 =====================
 
+* ENH: Add a Framewise Displacement calculation interface (https://github.com/nipy/nipype/pull/1604)
 * FIX: Use builtins open and unicode literals for py3 compatibility (https://github.com/nipy/nipype/pull/1572)
 * TST: reduce the size of docker images & use tags for images (https://github.com/nipy/nipype/pull/1564)
 * ENH: Implement missing inputs/outputs in FSL AvScale (https://github.com/nipy/nipype/pull/1563)

nipype/algorithms/confounds.py

@@ -0,0 +1,154 @@
+# -*- coding: utf-8 -*-
+# emacs: -*- mode: python; py-indent-offset: 4; indent-tabs-mode: nil -*-
+# vi: set ft=python sts=4 ts=4 sw=4 et:
+'''
+Algorithms to compute confounds in :abbr:`fMRI (functional MRI)`
+
+    Change directory to provide relative paths for doctests
+    >>> import os
+    >>> filepath = os.path.dirname(os.path.realpath(__file__))
+    >>> datadir = os.path.realpath(os.path.join(filepath, '../testing/data'))
+    >>> os.chdir(datadir)
+
+'''
+from __future__ import print_function, division, unicode_literals, absolute_import
+from builtins import str, zip, range, open
+
+import os.path as op
+import numpy as np
+
+from .. import logging
+from ..external.due import due, BibTeX
+from ..interfaces.base import (traits, TraitedSpec, BaseInterface,
+                               BaseInterfaceInputSpec, File, isdefined)
+IFLOG = logging.getLogger('interface')
+
+
+class FramewiseDisplacementInputSpec(BaseInterfaceInputSpec):
+    in_plots = File(exists=True, desc='motion parameters as written by FSL MCFLIRT')
+    radius = traits.Float(50, usedefault=True,
+                          desc='radius in mm to calculate angular FDs, 50mm is the '
+                               'default since it is used in Power et al. 2012')
+    out_file = File('fd_power_2012.txt', usedefault=True, desc='output file name')
+    out_figure = File('fd_power_2012.pdf', usedefault=True, desc='output figure name')
+    series_tr = traits.Float(desc='repetition time in sec.')
+    save_plot = traits.Bool(False, usedefault=True, desc='write FD plot')
+    normalize = traits.Bool(False, usedefault=True, desc='calculate FD in mm/s')
+    figdpi = traits.Int(100, usedefault=True, desc='output dpi for the FD plot')
+    figsize = traits.Tuple(traits.Float(11.7), traits.Float(2.3), usedefault=True,
+                           desc='output figure size')
+
+class FramewiseDisplacementOutputSpec(TraitedSpec):
+    out_file = File(desc='calculated FD per timestep')
+    out_figure = File(desc='output image file')
+    fd_average = traits.Float(desc='average FD')
+
+class FramewiseDisplacement(BaseInterface):
+    """
+    Calculate the :abbr:`FD (framewise displacement)` as in [Power2012]_.
+    This implementation reproduces the calculation in fsl_motion_outliers
+
+    .. [Power2012] Power et al., Spurious but systematic correlations in functional
+         connectivity MRI networks arise from subject motion, NeuroImage 59(3),
+         2012. doi:`10.1016/j.neuroimage.2011.10.018
+         <http://dx.doi.org/10.1016/j.neuroimage.2011.10.018>`_.
+
+
+    """
+
+    input_spec = FramewiseDisplacementInputSpec
+    output_spec = FramewiseDisplacementOutputSpec
+
+    references_ = [{
+        'entry': BibTeX("""\
+@article{power_spurious_2012,
+    title = {Spurious but systematic correlations in functional connectivity {MRI} networks \
+arise from subject motion},
+    volume = {59},
+    doi = {10.1016/j.neuroimage.2011.10.018},
+    number = {3},
+    urldate = {2016-08-16},
+    journal = {NeuroImage},
+    author = {Power, Jonathan D. and Barnes, Kelly A. and Snyder, Abraham Z. and Schlaggar, \
+Bradley L. and Petersen, Steven E.},
+    year = {2012},
+    pages = {2142--2154},
+}
+"""),
+        'tags': ['method']
+    }]
+
+    def _run_interface(self, runtime):
+        mpars = np.loadtxt(self.inputs.in_plots)  # mpars is N_t x 6
+        diff = mpars[:-1, :] - mpars[1:, :]
+        diff[:, :3] *= self.inputs.radius
+        fd_res = np.abs(diff).sum(axis=1)
+
+        self._results = {
+            'out_file': op.abspath(self.inputs.out_file),
+            'fd_average': float(fd_res.mean())
+        }
+        np.savetxt(self.inputs.out_file, fd_res)
+
+
+        if self.inputs.save_plot:
+            tr = None
+            if isdefined(self.inputs.series_tr):
+                tr = self.inputs.series_tr
+
+            if self.inputs.normalize and tr is None:
+                IFLOG.warn('FD plot cannot be normalized if TR is not set')
+
+            self._results['out_figure'] = op.abspath(self.inputs.out_figure)
+            fig = plot_confound(fd_res, self.inputs.figsize, 'FD', units='mm',
+                                series_tr=tr, normalize=self.inputs.normalize)
+            fig.savefig(self._results['out_figure'], dpi=float(self.inputs.figdpi),
+                        format=self.inputs.out_figure[-3:],
+                        bbox_inches='tight')
+            fig.clf()
+        return runtime
+
+    def _list_outputs(self):
+        return self._results
+
+
+def plot_confound(tseries, figsize, name, units=None,
+                  series_tr=None, normalize=False):
+    """
+    A helper function to plot :abbr:`fMRI (functional MRI)` confounds.
+    """
+    import matplotlib
+    matplotlib.use('Agg')
+    import matplotlib.pyplot as plt
+    from matplotlib.gridspec import GridSpec
+    from matplotlib.backends.backend_pdf import FigureCanvasPdf as FigureCanvas
+    import seaborn as sns
+
+    fig = plt.Figure(figsize=figsize)
+    FigureCanvas(fig)
+    grid = GridSpec(1, 2, width_ratios=[3, 1], wspace=0.025)
+    grid.update(hspace=1.0, right=0.95, left=0.1, bottom=0.2)
+
+    ax = fig.add_subplot(grid[0, :-1])
+    if normalize and series_tr is not None:
+        tseries /= series_tr
+
+    ax.plot(tseries)
+    ax.set_xlim((0, len(tseries)))
+    ylabel = name
+    if units is not None:
+        ylabel += (' speed [{}/s]' if normalize else ' [{}]').format(units)
+    ax.set_ylabel(ylabel)
+
+    xlabel = 'Frame #'
+    if series_tr is not None:
+        xlabel = 'Frame # ({} sec TR)'.format(series_tr)
+    ax.set_xlabel(xlabel)
+    ylim = ax.get_ylim()
+
+    ax = fig.add_subplot(grid[0, -1])
+    sns.distplot(tseries, vertical=True, ax=ax)
+    ax.set_xlabel('Frames')
+    ax.set_ylim(ylim)
+    ax.set_yticklabels([])
+    return fig

nipype/algorithms/misc.py

@@ -1157,6 +1157,10 @@ class SplitROIs(BaseInterface):
     """
     Splits a 3D image in small chunks to enable parallel processing.
     ROIs keep time series structure in 4D images.
+
+    Example
+    -------
+
     >>> from nipype.algorithms import misc
     >>> rois = misc.SplitROIs()
     >>> rois.inputs.in_file = 'diffusion.nii'

nipype/algorithms/tests/test_auto_FramewiseDisplacement.py

@@ -0,0 +1,43 @@
+# AUTO-GENERATED by tools/checkspecs.py - DO NOT EDIT
+from ...testing import assert_equal
+from ..confounds import FramewiseDisplacement
+
+
+def test_FramewiseDisplacement_inputs():
+    input_map = dict(figdpi=dict(usedefault=True,
+    ),
+    figsize=dict(usedefault=True,
+    ),
+    ignore_exception=dict(nohash=True,
+    usedefault=True,
+    ),
+    in_plots=dict(),
+    normalize=dict(usedefault=True,
+    ),
+    out_figure=dict(usedefault=True,
+    ),
+    out_file=dict(usedefault=True,
+    ),
+    radius=dict(usedefault=True,
+    ),
+    save_plot=dict(usedefault=True,
+    ),
+    series_tr=dict(),
+    )
+    inputs = FramewiseDisplacement.input_spec()
+
+    for key, metadata in list(input_map.items()):
+        for metakey, value in list(metadata.items()):
+            yield assert_equal, getattr(inputs.traits()[key], metakey), value
+
+
+def test_FramewiseDisplacement_outputs():
+    output_map = dict(fd_average=dict(),
+    out_figure=dict(),
+    out_file=dict(),
+    )
+    outputs = FramewiseDisplacement.output_spec()
+
+    for key, metadata in list(output_map.items()):
+        for metakey, value in list(metadata.items()):
+            yield assert_equal, getattr(outputs.traits()[key], metakey), value

nipype/algorithms/tests/test_confounds.py

@@ -0,0 +1,18 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+
+from nipype.testing import (assert_equal, example_data)
+from nipype.algorithms.confounds import FramewiseDisplacement
+import numpy as np
+from tempfile import mkdtemp
+from shutil import rmtree
+
+def test_fd():
+    tempdir = mkdtemp()
+    ground_truth = np.loadtxt(example_data('fsl_motion_outliers_fd.txt'))
+    fd = FramewiseDisplacement(in_plots=example_data('fsl_mcflirt_movpar.txt'),
+                               out_file=tempdir + '/fd.txt')
+    res = fd.run()
+    yield assert_equal, np.allclose(ground_truth, np.loadtxt(res.outputs.out_file)), True
+    yield assert_equal, np.abs(ground_truth.mean() - res.outputs.fd_average) < 1e-4, True
+    rmtree(tempdir)