Closed oesteban closed 1 year ago
Merging #166 (b97e55c) into master (b610a9a) will increase coverage by
0.01%
. The diff coverage is100.00%
.
@@ Coverage Diff @@
## master #166 +/- ##
==========================================
+ Coverage 98.60% 98.62% +0.01%
==========================================
Files 13 13
Lines 1217 1232 +15
Branches 184 187 +3
==========================================
+ Hits 1200 1215 +15
Misses 10 10
Partials 7 7
Flag | Coverage Δ | |
---|---|---|
travis | 96.91% <100.00%> (+0.12%) |
:arrow_up: |
unittests | 98.57% <100.00%> (+0.01%) |
:arrow_up: |
Flags with carried forward coverage won't be shown. Click here to find out more.
Impacted Files | Coverage Δ | |
---|---|---|
nitransforms/__init__.py | 100.00% <100.00%> (ø) |
|
nitransforms/manip.py | 100.00% <100.00%> (ø) |
|
nitransforms/nonlinear.py | 98.95% <100.00%> (+0.19%) |
:arrow_up: |
Continue to review full report at Codecov.
Legend - Click here to learn more
Δ = absolute <relative> (impact)
,ø = not affected
,? = missing data
Powered by Codecov. Last update b610a9a...b97e55c. Read the comment docs.
I'm confused, why a deformation field is not applicable? and conversely, why a displacements field is not composable?
ITK generates displacements fields and you can ask ANTS to generate "composite" transforms, which in effect are a new displacements field.
Similarly, if you use SPM, nonlinear transforms most often will be given in deformation convention (if I correctly understood J. Ashburner in our meeting in Boston back in 2018) and again, you can apply/compose as you wish.
Perhaps, the missing piece here is a convenience function to convert from deformation to displacements and define __matmul__
on the displacements class so that another displacements object is returned.
One more reason to use the proposed design is #155. With the current implementation, that one becomes trivial for both types of fields.
My understanding here is that DeformationFieldTransform
has an identity during composition and application, while DisplacementsFieldTransform
can have an identity during application (and it's essentially the same), but that identity does not hold during composition.
I have not tested the following, but here's my rough expectation from reading this code:
>>> dfm = DeformationFieldTransform()
>>> disp = DisplacementsFieldTransform()
Application
>>> x, y, z = some_coord
>>> np.array_equal(dfm.map([x, y, z]), [x, y, z]))
True
>>> np.array_equal(disp.map([x, y, z]), [x, y, z]))
True
Composition
>>> dfm @ dfm == dfm
True
>>> disp @ disp == disp
False
Maybe a better way to put this is that a DeformationFieldTransform
seems like the correct internal representation, and that we should not have a DisplacementsFieldTransform
except possibly as an implementation detail for certain formats.
Okay, I think I now understand what you mean.
If you assume that:
>>> dfm = DeformationFieldTransform()
>>> disp = DisplacementsFieldTransform()
will create identity transforms for both structures then:
dfm._field
contains an $(S_x, S_y, S_z, D)$ array where each voxel contains the three coordinates of that voxel in physical space.disp._displacements
is all zeros, and disp._field == dfm._field
.Then, mapping is equivalent for both.
Regarding composition, in both cases:
>>> dfm @ dfm == dfm
True
>>> disp @ disp == disp
True
iff they are identity.
For both of them, if they are not identity, then both satisfy:
>>> dfm @ dfm == dfm
False
>>> disp @ disp == disp
False
which will become more apparent with the implementation of #155. In the case of disp
I think you were already hinting at this, but it is exactly the same for dfm
.
Perhaps you would feel more comfortable with a single representation structure:
class DenseFieldTransform(TransformBase):
def __init__(self, deformation=None, displacements=None, reference=None):
...
where:
self._field
) or displacements (self._displacements
) are calculated from the provided argumentNone
, an identity field is generated.Does this make sense?
Okay, I sketched out an alternative API. Let me know what you think. Still, I need to implement the identity, but should be easy to do.
Okay, reading the new API, I think some things that were unclear are clicking. Now I'm just left to wonder: Why do we want to have the _displacements
pre-calculated at all? Why do we not just use:
class DenseFieldTransform(TransformBase):
def __init__(self, deltas, reference=None):
self._deltas = deltas
self._reference = reference
def __matmul__(self, other):
deltas = b.map(
self._deltas.reshape((-1, self._deltas.shape[-1]))
).reshape(self._deltas.shape)
return DenseFieldTransform(deltas, reference=self.reference)
def map(self, coords, inverse=False):
indices = calculate_indices(coords)
return coords + self._deltas[indices]
@classmethod
def from_deformation_field(cls, deformations, reference=None):
if reference is None:
reference = ...
deltas = deformations - reference.ndcoords.T.reshape(deformations.shape)
return cls(deltas, reference=reference)
Finally, are we sure that we want to use __matmul__
and not __rmatmul__
? We would need to think of a transformation that would be both simple and obviously not commutative to verify that this goes in the direction we expect, but I can't immediately think of one.
Apologies if I'm jumping all over the place. Trying to understand this somewhat quickly and may be flailing a bit...
Why do we want to have the
_displacements
pre-calculated at all?
I think this is just for reproducibility purposes -- to be able to write it out to disk (possibly in a different format such as X5) without alteration of the array contents.
Why do we not just use:
As commented in the code, AFAIK the deformations are numerically more stable and that's why I made them the default.
Finally, are we sure that we want to use
__matmul__
and not__rmatmul__
? We would need to think of a transformation that would be both simple and obviously not commutative to verify that this goes in the direction we expect, but I can't immediately think of one.
I opened an issue to remind me of this good point.
Apologies if I'm jumping all over the place. Trying to understand this somewhat quickly and may be flailing a bit...
On the contrary, this is very helpful. Thanks for your time.
I'll use your _deltas
convention because I like it better than _displacements
. Other than that, I will call this PR finished and move on. Thanks a lot for the feedback. Happy to come back to any comments anytime.
This PR lays the ground for future work on #56, and #89, by defining the matrix multiplication operator on field-based transforms.