Affine registration PR1: Transforms. #574

omarocegueda · 2015-02-05T23:41:19Z

Based on some performance comparisons between ANTS2, Flirt and Nipy for affine registration with mutual information, we have observed that ANTS2 consistently performs better in terms of registration accuracy. The following graph depicts the Jaccard index of 31 manually annotated anatomical regions averaged over 306 registrations (1=perfect overlap, 0=no overlap):
http://www.cimat.mx/~omar/imgshare/mi32_flirt_ants2_nipy.png

This motivated us to implement the same registration strategy in Cython, so we can easily use it with our non-linear registration module and at some point modify the algorithms/metrics more easily for DWMRI. This is the current performance of this Cython implementation:
http://www.cimat.mx/~omar/imgshare/ANTS_vs_Dipy.png
I am investigating the cause of the small differences but at some point we should reach the same performance.

In an effort to make the review process (a bit) easier, this part of the code only contains the computation of the affine transforms and their Jacobians.

stefanv · 2015-02-06T02:12:15Z

dipy/align/transforms.pyx

+    return 0
+
+
+cdef void _dot_prod(double[:,:] A, double[:,:] B, double[:,:] C):


Not tested or used?

stefanv · 2015-02-06T02:13:22Z

@omarocegueda I like your structured approach to measuring these results. Can we include some form of the performance (in terms of accuracy) tests in the code base?

omarocegueda · 2015-02-06T17:49:48Z

Thank you @stefanv!,
sure, it is actually issue #441, the problem is that the full evaluation may take a long time (especially for non-linear registration). The way I do it right now is to create one job for each registration and send the (306) jobs to a cluster. Maybe we can start by adding only some evaluation tools (dipy.align.eval?) to the align module (code to compute the Jaccard indices, and creating the graphs), and assume that the users somehow ran all the registrations they need to evaluate.

matthew-brett · 2015-02-19T22:54:13Z

dipy/align/tests/test_transforms.py

@@ -0,0 +1,233 @@
+from dipy.align.transforms import *


Please import specific names (best) or import the module (import dipy.align.transforms as dit) and use with dit. namespace. The * makes it difficult to spot bad imports, among other things.

matthew-brett · 2015-02-19T23:31:18Z

Sorry - this might be silly - but how about having a set of strategy objects like:

cpdef class RegMetric:
    cdef get_jacobian(double[:], double[:], double[:,:]) nogil:
    cdef param_to_matrix_function(double[:], double[:]) nogil

Instantiate all the ones you need for 2D and 3D, different transformations e.g:

regmetrics = {}
regmetrics[(TRANSLATION, 2)] = Translation2D()

etc, and then do regmetrics[(ttype, dims)] to get the strategy you need?

… of docstring in test functions. Document return value of Jacobian functions. Remove 'with nogil' clause when calling a nogil function. Associated to --> Associated with. Removed unused '_dot_prod' function.

…d argument handling.

matthew-brett · 2015-02-26T01:21:52Z

dipy/align/transforms.pyx

+        J : array, shape (dim, n)
+            Jacobian matrix of the transform with parameters theta at point x
+        """
+        n = len(theta)


Python call? Maybe use theta.shape[0]

Garyfallidis · 2015-02-26T02:11:43Z

Hi @omarocegueda, I see only 3 files in this PR if this is going to be only about the transforms, you should change the title of the PR. Travis is failing also. Can you correct that?

…tead of iteritems (deprecated in Python 3)

omarocegueda · 2015-02-26T03:39:38Z

Hi @Garyfallidis!, yes, this PR only contains the Transforms part. I just changed the title. Travis was failing because "iteritems" was deprecated in Python3, now using "items".

Garyfallidis · 2015-02-26T18:28:31Z

Do you still need the Work In Progress (WIP) flag on the title? When we see WIP it means that we need to wait until this PR is ready for a review. You can change WIP to MRG if you think this is ready to be merged on your side or just remove the WIP.

Garyfallidis · 2015-02-26T20:02:15Z

dipy/align/transforms.pyx

+        theta : array, shape(2,)
+            buffer to write the parameters of the 2D translation transform
+        """
+        theta[:2] = 0


Is this really applied properly in Cython? I mean in a cdef with nogil.
So, this will do theta[0] = 0 and theta[1] = 1? No for loop needed?

I think so : http://docs.cython.org/src/userguide/memoryviews.html

# NumPy-style syntax for assigning a single value to all elements. narr_view[:, :, :] = 3

Cool! I forgot about that. Thank you @matthew-brett.

matthew-brett · 2015-02-26T20:46:12Z

Omar - do you have a feel for how good the test coverage is for the code?

Garyfallidis · 2015-02-26T20:52:48Z

dipy/align/transforms.pyx

+                  ( cc * sb + sc * sa * cb ) * pz
+        J[2, 2] = 0
+
+        J[0,3:6] = 0


space missing J[0, 3:6]

Garyfallidis · 2015-02-26T22:22:05Z

Hey @omarocegueda please address mine and Matthew's comments. Apart from that. This PR looks very nice and I believe is very close to be merged. Thank you for all the good work and for separating your larger PR in smaller and easier to follow PRs.

…ions

omarocegueda · 2015-02-27T00:05:05Z

Hi @matthew-brett!, I tried to use your find-cpdef-tests, but I am having some weird errors, I need to see what went wrong. However, we have close to full coverage: in the tests I iterate over all transforms calling all of their methods, so all methods are called at least once. I also check the ValueError exceptions. I think the only missing case is what happens when the user attempts to use the base class Transform instead of one of its derived classes (which is undefined), what do you think would be a good way to handle that case?(an abstract class would be ideal, so it couldn't be instantiated, afaik there is no such a thing as an "abstract struct"). Right now it will raise ValueError when attempting to get its Jacobian or its identity parameters or its matrix representation, and it returns -1 when queried about its dimension and number of parameters.

matthew-brett · 2015-02-27T00:07:55Z

Yes, my function was only designed for Cython functions, so won't be much use in this case, with lots of classes.

It sounds as if the error modes for the Transform are OK - maybe a few tests for those just to remind ourselves that there are errors and the reason for them.

omarocegueda · 2015-02-27T00:30:19Z

Alright, we should have full coverage now.
Thanks!

matthew-brett · 2015-02-27T00:31:19Z

Good to go from my point of view - thanks Omar.

Garyfallidis · 2015-02-27T18:43:28Z

Thx Omar. Let's move to the second PR.

Affine registration PR1: Transforms.

Add code for affine transformations and their Jacobians.

170e123

stefanv reviewed Feb 6, 2015
View reviewed changes

dipy/align/transforms.pyx

return 0

cdef void _dot_prod(double[:,:] A, double[:,:] B, double[:,:] C):

Copy link

Contributor

stefanv Feb 6, 2015

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Not tested or used?

omarocegueda changed the title ~~WIP: mlti-resolution affine registration with mutual information.~~ WIP: multi-resolution affine registration with mutual information. Feb 6, 2015

matthew-brett reviewed Feb 19, 2015
View reviewed changes

Import individual symbols from transforms module. Use comment instead…

937748f

… of docstring in test functions. Document return value of Jacobian functions. Remove 'with nogil' clause when calling a nogil function. Associated to --> Associated with. Removed unused '_dot_prod' function.

omarocegueda force-pushed the imaffine_pr1 branch from 14225b6 to 937748f Compare February 23, 2015 17:50

Validate arguments, raise ValueError if invalid. Add tests for invali…

a8330ea

…d argument handling.

matthew-brett reviewed Feb 26, 2015
View reviewed changes

Using .shape[0] instead of len(.) to avoid python call. Use items ins…

10152af

…tead of iteritems (deprecated in Python 3)

omarocegueda changed the title ~~WIP: multi-resolution affine registration with mutual information.~~ WIP: Affine registration PR1: Transforms. Feb 26, 2015

Garyfallidis reviewed Feb 26, 2015
View reviewed changes

dipy/align/transforms.pyx

( cc * sb + sc * sa * cb ) * pz

J[2, 2] = 0

J[0,3:6] = 0

Copy link

Contributor

Garyfallidis Feb 26, 2015

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

space missing J[0, 3:6]

omarocegueda changed the title ~~WIP: Affine registration PR1: Transforms.~~ Affine registration PR1: Transforms. Feb 26, 2015

Use np.random.RandomState instead of np.random.rand, and PEP8 correct…

7780dff

…ions

omarocegueda force-pushed the imaffine_pr1 branch from 3902149 to 7780dff Compare February 26, 2015 22:44

omarocegueda force-pushed the imaffine_pr1 branch from 8e7e8e4 to b217bac Compare February 27, 2015 00:29

Add tests for invalid use of Transform

5f87d9e

omarocegueda force-pushed the imaffine_pr1 branch from b217bac to 5f87d9e Compare February 27, 2015 00:32

Garyfallidis added a commit that referenced this pull request Feb 27, 2015

Merge pull request #574 from omarocegueda/imaffine_pr1

7e3b95e

Affine registration PR1: Transforms.

Garyfallidis merged commit 7e3b95e into dipy:master Feb 27, 2015

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Affine registration PR1: Transforms. #574

Affine registration PR1: Transforms. #574

omarocegueda commented Feb 5, 2015

stefanv Feb 6, 2015

stefanv commented Feb 6, 2015

omarocegueda commented Feb 6, 2015

matthew-brett Feb 19, 2015

matthew-brett commented Feb 19, 2015

matthew-brett Feb 26, 2015

Garyfallidis commented Feb 26, 2015

omarocegueda commented Feb 26, 2015

Garyfallidis commented Feb 26, 2015

Garyfallidis Feb 26, 2015

matthew-brett Feb 26, 2015

Garyfallidis Feb 26, 2015

matthew-brett commented Feb 26, 2015

Garyfallidis Feb 26, 2015

Garyfallidis commented Feb 26, 2015

omarocegueda commented Feb 27, 2015

matthew-brett commented Feb 27, 2015

omarocegueda commented Feb 27, 2015

matthew-brett commented Feb 27, 2015

Garyfallidis commented Feb 27, 2015

		return 0


		cdef void _dot_prod(double[:,:] A, double[:,:] B, double[:,:] C):

Affine registration PR1: Transforms. #574

Affine registration PR1: Transforms. #574

Conversation

omarocegueda commented Feb 5, 2015

stefanv Feb 6, 2015

Choose a reason for hiding this comment

stefanv commented Feb 6, 2015

omarocegueda commented Feb 6, 2015

matthew-brett Feb 19, 2015

Choose a reason for hiding this comment

matthew-brett commented Feb 19, 2015

matthew-brett Feb 26, 2015

Choose a reason for hiding this comment

Garyfallidis commented Feb 26, 2015

omarocegueda commented Feb 26, 2015

Garyfallidis commented Feb 26, 2015

Garyfallidis Feb 26, 2015

Choose a reason for hiding this comment

matthew-brett Feb 26, 2015

Choose a reason for hiding this comment

Garyfallidis Feb 26, 2015

Choose a reason for hiding this comment

matthew-brett commented Feb 26, 2015

Garyfallidis Feb 26, 2015

Choose a reason for hiding this comment

Garyfallidis commented Feb 26, 2015

omarocegueda commented Feb 27, 2015

matthew-brett commented Feb 27, 2015

omarocegueda commented Feb 27, 2015

matthew-brett commented Feb 27, 2015

Garyfallidis commented Feb 27, 2015