Merge pull request #51 from oesteban/enh/phdiff-better-demean-demode

 ENH: Do not use legacy demean function from old nipype workflows

.circleci/config.yml

@@ -224,20 +224,12 @@ jobs:
           paths:
             - /tmp/docker
       - run:
-          name: Set-up a Docker registry
-          command: |
-            docker run -d -p 5000:5000 --restart=always --name=registry \
-                -v /tmp/docker:/var/lib/registry registry:2
-      - run:
-          name: Pull images from local registry
-          command: |
-            docker pull localhost:5000/sdcflows
-            docker tag localhost:5000/sdcflows poldracklab/sdcflows:latest
-      - run:
-          name: Refresh work directory?
+          name: Refreshing cached intermediate results
           command: |
+            cd /tmp/src/sdcflows
+            COMMIT_MSG=$( git log --format=oneline -n 1 $CIRCLE_SHA1 )
             set +e
-            do_refresh="$( git log --format=oneline -n 1 $CIRCLE_SHA1 | grep -i -E '\[fresh[ _]?workdir\]' )"
+            do_refresh="$( echo "${COMMIT_MSG}" | grep -i -E '\[fresh[ _]?workdir\]' )"
             set -e
             if [[ "x${do_refresh}" = "x" ]]; then
               echo "Did not refresh the workdir."
@@ -247,6 +239,23 @@ jobs:
               cd /tmp/work
               tar xzfv /tmp/data/workdir.tar.gz
             fi
+
+            wipe_dir=$( echo "${COMMIT_MSG}" | sed -n 's/.*\[wipe \([a-zA-Z0-9_\*]*\)\].*/\1/p' )
+            if [[ "x${wipe_dir}" != "x" ]]; then
+              path=/tmp/work/${wipe_dir}
+              echo "Found tag [wipe ${wipe_dir}] - clearing up $path ..."
+              rm -rf ${path}
+            fi
+      - run:
+          name: Set-up a Docker registry
+          command: |
+            docker run -d -p 5000:5000 --restart=always --name=registry \
+                -v /tmp/docker:/var/lib/registry registry:2
+      - run:
+          name: Pull images from local registry
+          command: |
+            docker pull localhost:5000/sdcflows
+            docker tag localhost:5000/sdcflows poldracklab/sdcflows:latest
       - run:
           name: Run tests
           no_output_timeout: 2h

sdcflows/workflows/phdiff.py

@@ -19,8 +19,7 @@
 
 from nipype.interfaces import ants, fsl, utility as niu
 from nipype.pipeline import engine as pe
-from niflow.nipype1.workflows.dmri.fsl.utils import (
-    demean_image, cleanup_edge_pipeline)
+from niflow.nipype1.workflows.dmri.fsl.utils import cleanup_edge_pipeline
 from niworkflows.engine.workflows import LiterateWorkflow as Workflow
 from niworkflows.interfaces.images import IntraModalMerge
 from niworkflows.interfaces.masks import BETRPT
@@ -125,7 +124,7 @@ def init_phdiff_wf(omp_nthreads, name='phdiff_wf'):
     denoise = pe.Node(fsl.SpatialFilter(operation='median', kernel_shape='sphere',
                                         kernel_size=3), name='denoise')
 
-    demean = pe.Node(niu.Function(function=demean_image), name='demean')
+    demean = pe.Node(niu.Function(function=_demean), name='demean')
 
     cleanup_wf = cleanup_edge_pipeline(name="cleanup_wf")
 
@@ -171,3 +170,38 @@ def _recenter(in_file):
                                newpath=getcwd())
     nb.Nifti1Image(data, nii.affine, nii.header).to_filename(out_file)
     return out_file
+
+
+def _demean(in_file, in_mask=None, usemode=True):
+    """
+    Subtract the median (since it is robuster than the mean) from a map.
+
+    Parameters
+    ----------
+    usemode : bool
+        Use the mode instead of the median (should be even more robust
+        against outliers).
+
+    """
+    from os import getcwd
+    import numpy as np
+    import nibabel as nb
+    from nipype.utils.filemanip import fname_presuffix
+
+    nii = nb.load(in_file)
+    data = nii.get_fdata(dtype='float32')
+
+    msk = np.ones_like(data, dtype=bool)
+    if in_mask is not None:
+        msk[nb.load(in_mask).get_fdata(dtype='float32') < 1e-4] = False
+
+    if usemode:
+        from scipy.stats import mode
+        data[msk] -= mode(data[msk], axis=None)[0][0]
+    else:
+        data[msk] -= np.median(data[msk], axis=None)
+
+    out_file = fname_presuffix(in_file, suffix='_demean',
+                               newpath=getcwd())
+    nb.Nifti1Image(data, nii.affine, nii.header).to_filename(out_file)
+    return out_file