Merge branch 'fix-shoshber-compcor1594' of https://github.com/oesteban/nipype into oesteban-fix-shoshber-compcor1594

Author: Shoshana Berleant

Makefile

@@ -53,7 +53,10 @@ inplace:
 	$(PYTHON) setup.py build_ext -i
 
 test-code: in
-	$(NOSETESTS) -s nipype --with-doctest
+#	$(NOSETESTS) -s nipype/algorithms/compcor.py --with-doctest -v
+	$(NOSETESTS) -s nipype/algorithms/tests/test_tsnr.py #--with-coverage
+	$(NOSETESTS) -s nipype/algorithms/tests/test_compcor.py #--with-coverage
+#	$(NOSETESTS) -s nipype/algorithms/misc.py --with-doctest -v
 
 test-doc:
 	$(NOSETESTS) -s --with-doctest --doctest-tests --doctest-extension=rst \

nipype/algorithms/compcor.py

@@ -123,6 +123,7 @@ def _run_interface(self, runtime):
         mask_file = 'mask.nii'
         nb.nifti1.save(nb.Nifti1Image(mask, np.eye(4)), mask_file)
         self.inputs.mask_file = 'mask.nii'
+
         super(TCompCor, self)._run_interface(runtime)
         return runtime
 

nipype/algorithms/tsnr.py

@@ -60,8 +60,11 @@ def _run_interface(self, runtime):
         img = nb.load(self.inputs.in_file[0])
         header = img.header.copy()
         vollist = [nb.load(filename) for filename in self.inputs.in_file]
-        data = np.concatenate([vol.get_data().reshape(
-            vol.get_shape()[:3] + (-1,)) for vol in vollist], axis=3)
+        data = np.concatenate(
+            [vol.get_data().reshape(
+                vol.get_shape()[:3] + (-1,))
+             for vol in vollist],
+                              axis=3)
         data = np.nan_to_num(data)
 
         if data.dtype.kind == 'i':
@@ -98,17 +101,25 @@ def regress_poly(degree, data):
     ''' returns data with degree polynomial regressed out.
     The last dimension (i.e. data.shape[-1]) should be time.
     '''
-    timepoints = data.shape[-1]
+    datashape = data.shape
+    timepoints = datashape[-1]
+
+    # Rearrange all voxel-wise time-series in rows
+    data = data.reshape((-1, timepoints))
+
+    # Generate design matrix
     X = np.ones((timepoints, 1))
     for i in range(degree):
         polynomial_func = legendre(i+1)
         value_array = np.linspace(-1, 1, timepoints)
         X = np.hstack((X, polynomial_func(value_array)[:, newaxis]))
 
-    betas = np.dot(np.linalg.pinv(X), np.rollaxis(data, 3, 2))
-    datahat = np.rollaxis(np.dot(X[:, 1:],
-                                 np.rollaxis(
-                                     betas[1:, :, :, :], 0, 3)),
-                          0, 4)
+    # Calculate coefficients
+    betas = np.linalg.pinv(X).dot(data.T)
+
+    # Estimation
+    datahat = X[:, 1:].dot(betas[1:, ...]).T
     regressed_data = data - datahat
-    return regressed_data
+
+    # Back to original shape
+    return regressed_data.reshape(datashape)