ENH: A script which serves as a wrapper to FmriRealign4d #265

Merged
merged 21 commits into from Apr 15, 2014

Projects

None yet

5 participants

@arokem
Member
arokem commented Mar 18, 2013

This is a little script to do argument parsing, file-handling, figure generation, around the 4D motion correction/realignment functionality. Emulates the input/output of FSL's mcflirt and also generates a figure with the parameters.

@arokem
Member
arokem commented Mar 18, 2013

I should say - this is still WIP - please don't merge yet. Thought I would put it up here to give you all a heads up that I am working on this.

@alexis-roche
Member

Thank you Ariel for doing this, it's actually very useful. Your code looks nice to me, but I would like to point you to the slice_order thread I just started on the nipy list - which most likely has no consequence for this PR but means the underlying FmriRealign4D code needs changes.

@arokem
Member
arokem commented Apr 14, 2013

OK - this is now ready to go, I think. I've added testing on the test data that is in nipy.testing and I've also tried it out on some of my own data and it seems to be doing the right thing. If others could give it a whirl, that would be great, of course.

@alexis-roche
Member

Please beware that, as of #266, the slice_order argument of FmriRealign4D is deprecated and the preferred way to declare slice acquisition timing is via a slice_times argument.

@arokem
Member
arokem commented Apr 15, 2013

Thanks @alexis-roche - I will keep track of what ends up happening there and refactor here as necessary.

@matthew-brett
Member

Hi Ariel - very sorry to be slow to get to this. There's a lot of code in the script. Can you move any code that is not script-specific out into the library part of the code?

@arokem
Member
arokem commented Apr 26, 2013

No worries - I am also moving along rather slowly with this. Thanks for
taking a look.

I will refactor to put stuff in the library. Any ideas where you think I
should put these things?

On Fri, Apr 26, 2013 at 9:34 AM, Matthew Brett notifications@github.comwrote:

Hi Ariel - very sorry to be slow to get to this. There's a lot of code in
the script. Can you move any code that is not script-specific out into the
library part of the code?


Reply to this email directly or view it on GitHubhttps://github.com/nipy/nipy/pull/265#issuecomment-17085302
.

@satra
Member
satra commented Apr 26, 2013

just a note over here - the angles from this algorithm are not euclidean angles and in a script such as this it might be useful to point out.

@arokem
Member
arokem commented Apr 26, 2013

Oh - I am not sure that I know what that means. Does this refer to the fact
that the angles are given in radians, rather than angles, or is there
something else that I am unaware of?

On Fri, Apr 26, 2013 at 11:31 AM, Satrajit Ghosh
notifications@github.comwrote:

just a note over here - the angles from this algorithm are not euclidean
angles and in a script such as this it might be useful to point out.


Reply to this email directly or view it on GitHubhttps://github.com/nipy/nipy/pull/265#issuecomment-17091993
.

@matthew-brett
Member

For location - how about nipy/algorithms/registration/scripting.py ?

@alexis-roche
Member

@satra: you mean Euler angles? You're right, it would be worth mentioning that the 3 rotation parameters output by Realign4D are not Euler angles (as in SPM and other common packages). They instead correspond to a rotation vector, that is, the rotation angle times the rotation axis normalized to unit norm:

http://en.wikipedia.org/wiki/Axis-angle_representation

It's just another representation of a rotation, neither better nor worse than Euler angles.

@satra
Member
satra commented Apr 26, 2013

thanks @alexis-roche - indeed that's i meant! :)

@matthew-brett
Member

Ariel - sorry - would you consider refactoring this against the current interface? That would be great help.

@arokem
Member
arokem commented Mar 26, 2014

Yeah - I'll get on it.

On Wed, Mar 26, 2014 at 10:40 AM, Matthew Brett notifications@github.comwrote:

Ariel - sorry - would you consider refactoring this against the current
interface? That would be great help.


Reply to this email directly or view it on GitHubhttps://github.com/nipy/nipy/pull/265#issuecomment-38715256
.

@arokem
Member
arokem commented Mar 26, 2014

Following @matthew-brett's recent (well, less than a year) comment, I have tried moving stuff into nipy/algorithms/registration/scripting.py in that last commit, but I can't import nipy.algorithms.registration.scripting. Any idea what I'm doing wrong?

@arokem
Member
arokem commented Mar 27, 2014

OK - I think I got it now. Let's see what Travis says.

@arokem
Member
arokem commented Mar 27, 2014

Oh dear. MPL's not a requirement, huh?

@matthew-brett
Member

No - did you see e.g nipy/tests/test_scripts.py? There's some stuff there
for dealing with MPL as an optional dependency.

On Wed, Mar 26, 2014 at 5:49 PM, Ariel Rokem notifications@github.comwrote:

Oh dear. MPL's not a requirement, huh?

Reply to this email directly or view it on GitHubhttps://github.com/nipy/nipy/pull/265#issuecomment-38758226
.

@arokem
Member
arokem commented Mar 27, 2014

This seems to have done the trick.

@satra
Member
satra commented Mar 27, 2014

looks good. one question i have is about adding the last value from a previous run to the transforms. aren't the transforms already relative to a "mean" position across runs?

also it might be useful to point out what the different parameters are being written out in the motion parameter file. (for example, to translate those parameters to a matrix folks will have to use nipy's to_matrix44 function - the spm or fsl like transforms won't work.)

@arokem
Member
arokem commented Mar 28, 2014

Thanks for the comments!

On Wed, Mar 26, 2014 at 8:18 PM, Satrajit Ghosh notifications@github.comwrote:

looks good. one question i have is about adding the last value from a
previous run to the transforms. aren't the transforms already relative to a
"mean" position across runs?

Yep - good catch. For some reason, I was using the
_within_runs_transforms attribute, so all of that acrobatics was needed,
but using the _transforms instead simplifies things substantially.
Comment: they all seem to be aligned to the first image in the first scan,
which is the desired behavior from my standpoint.

also it might be useful to point out what the different parameters are
being written out in the motion parameter file. (for example, to translate
those parameters to a matrix folks will have to use nipy's to_matrix44function - the spm or fsl like transforms won't work.)

I am not exactly sure what you mean, but the transforms are now returned
from this function and people can then do what they want with them (e.g.
use their as_affine attributes). I am guessing this is helpful?

Reply to this email directly or view it on GitHubhttps://github.com/nipy/nipy/pull/265#issuecomment-38765377
.

@coveralls

Coverage Status

Coverage remained the same when pulling affea46 on arokem:4drealign_wrapper into daaf644 on nipy:master.

@satra
Member
satra commented Mar 31, 2014

also it might be useful to point out what the different parameters are being written out in the motion parameter file. (for example, to translate those parameters to a matrix folks will have to use nipy's to_matrix44 function - the spm or fsl like transforms won't work.)

I am not exactly sure what you mean, but the transforms are now returned
from this function and people can then do what they want with them (e.g.
use their as_affine attributes). I am guessing this is helpful?

right - i was mainly pointing the executable script. that script will write out a file called xxx.par and most people will assume these are translations and rotations and apply the standard transforms to get at the affine matrix. but that won't work because the angles are not Euler angles. they would have to use nipy's to_matrix44. does this make sense?

@satra
Member
satra commented Mar 31, 2014

oh let me take that back - on closer inspection it doesn't save the translation and rotation, simply plots them. it might be good to save the transformation matrices then or the parameters, because they can then be used in later modeling stages.

@arokem
Member
arokem commented Mar 31, 2014

How's this?

On Sun, Mar 30, 2014 at 6:25 PM, Satrajit Ghosh notifications@github.comwrote:

oh let me take that back - on closer inspection it doesn't save the
translation and rotation, simply plots them. it might be good to save the
transformation matrices then or the parameters, because they can then be
used in later modeling stages.


Reply to this email directly or view it on GitHubhttps://github.com/nipy/nipy/pull/265#issuecomment-39047081
.

@satra satra and 1 other commented on an outdated diff Mar 31, 2014
scripts/nipy_4d_realign
+
+parser.add_argument('--slice_dim', type=int, metavar='Int', help="""Integer
+denoting the axis in `images` that is the slice axis. In a 4D image, this will
+often be axis = 2 (default).""", default=2)
+
+parser.add_argument('--slice_dir', type=int, metavar='Int', help=""" 1 if the
+slices were acquired slice 0 first (default), slice -1 last, or -1 if acquire slice -1 first, slice 0 last.""", default=1)
+
+parser.add_argument('--apply', type=bool, metavar='Bool',
+ help="""Whether to apply the realignment and save output files (default), or just estimate the realignment parameters and save those in '.par' files. {True, False}. Default: True""", default=True)
+
+parser.add_argument('--make_figure', type=bool, metavar='Bool',
+ help="""Whether to generate a '.png' figure with the motion parameters across runs. {True, False}. Default: False """, default=False)
+
+parser.add_argument('--save_params', type=bool, metavar='Bool',
+ help="""Whether to save the motion corrections parameters (3 rotations, 3 translations). {True, False}. Default: False """, default=False)
@satra
satra Mar 31, 2014 Member

Perhaps this help should note:
These rotations are not Euler angles, but a rotation vector. Please use nipyto_matrix44to convert to affine matrix.

@arokem
arokem Mar 31, 2014 Member

Yeah - I'm happy to add that.

On Mon, Mar 31, 2014 at 2:45 PM, Satrajit Ghosh notifications@github.comwrote:

In scripts/nipy_4d_realign:

+parser.add_argument('--slice_dim', type=int, metavar='Int', help="""Integer
+denoting the axis in images that is the slice axis. In a 4D image, this will
+often be axis = 2 (default).""", default=2)
+
+parser.add_argument('--slice_dir', type=int, metavar='Int', help=""" 1 if the
+slices were acquired slice 0 first (default), slice -1 last, or -1 if acquire slice -1 first, slice 0 last.""", default=1)
+
+parser.add_argument('--apply', type=bool, metavar='Bool',

  •            help="""Whether to apply the realignment and save output files (default), or just estimate the realignment parameters and save those in '.par' files. {True, False}. Default: True""", default=True)
    
    +parser.add_argument('--make_figure', type=bool, metavar='Bool',
  •            help="""Whether to generate a '.png' figure with the motion parameters across runs. {True, False}. Default: False """, default=False)
    
    +parser.add_argument('--save_params', type=bool, metavar='Bool',
  •            help="""Whether to save the motion corrections parameters (3 rotations, 3 translations). {True, False}. Default: False """, default=False)
    

Perhaps this help should note:
These rotations are not Euler angles, but a rotation vector. Please use
nipyto_matrix44to convert to affine matrix.


Reply to this email directly or view it on GitHubhttps://github.com/nipy/nipy/pull/265/files#r11138616
.

@matthew-brett
Member

Maybe it would be good to output Euler angles, for the sake of compatibility?

@arokem
Member
arokem commented Apr 1, 2014

Compatibility with what? I'm sorry, but I don't know how I would do this
conversion - could anyone give me a hint?

On Mon, Mar 31, 2014 at 6:44 PM, Matthew Brett notifications@github.comwrote:

Maybe it would be good to output Euler angles, for the sake of
compatibility?


Reply to this email directly or view it on GitHubhttps://github.com/nipy/nipy/pull/265#issuecomment-39162631
.

@matthew-brett
Member

Compatibility with SPM at least - I don't know what FSL plots, but I bet it's Euler angles too.

If you can get the 4x4 affine, then you can get the Euler angles with something like:

import numpy as np
import numpy.linalg as npl

from nibabel.eulerangles import mat2euler, euler2mat

from numpy.testing import assert_almost_equal

def aff2euler(affine):
    """ Return Euler angles from 4 x 4 `affine`
    """
    return mat2euler(aff2rot_zooms(affine)[0])


def aff2rot_zooms(affine):
    RZS = affine[:3, :3]
    zooms = np.sqrt(np.sum(RZS * RZS, axis=0))
    RS = RZS / zooms
    # Adjust zooms to make RS correspond (below) to a true
    # rotation matrix.
    if npl.det(RS) < 0:
        zooms[0] *= -1
        RS[:,0] *= -1
    # retrieve rotation matrix from RS with polar decomposition.
    # Discard shears
    P, S, Qs = npl.svd(RS)
    R = np.dot(P, Qs)
    return R, zooms


def test_aff2euler():
    xr = 0.1
    yr = -1.3
    zr = 3.1
    scales = (2.1, 3.2, 4.4)
    R = euler2mat(xr, yr, zr).dot(np.diag(scales))
    aff = np.eye(4)
    aff[:3, :3] = R
    aff[:3, 3] = [11, 12, 13]
    assert_almost_equal(aff2euler(aff), (xr, yr, zr))
@matthew-brett
Member

Ariel - where are we here?

@matthew-brett

zooms is a length 3 vector with voxel sizes

@arokem
Member
arokem commented Apr 10, 2014

We are trying to understand why this test fails:

https://travis-ci.org/nipy/nipy/jobs/22543936

The command-line input to make_figure is False:

https://github.com/arokem/nipy/blob/4drealign_wrapper/nipy/tests/test_scripts.py#L179

So - it shouldn't be going into that error statement. I don't understand it
and am still trying to reproduce this on my machine. I have about -60
minutes/day to work on this, considering how important this is to my
current research.

On Wed, Apr 9, 2014 at 7:12 PM, Matthew Brett notifications@github.comwrote:

Ariel - where are we here?


Reply to this email directly or view it on GitHubhttps://github.com/nipy/nipy/pull/265#issuecomment-40036607
.

@arokem
Member
arokem commented Apr 14, 2014

Alright - this is finally passing the tests. Is it ready to go, or does anyone have any more comments?

@matthew-brett
Member

All OK with me - any more comments before I merge?

@satra
Member
satra commented Apr 14, 2014

+1 for MRG

@matthew-brett matthew-brett merged commit f2fd307 into nipy:master Apr 15, 2014

1 check passed

continuous-integration/travis-ci The Travis CI build passed
Details
@arokem arokem deleted the arokem:4drealign_wrapper branch Apr 15, 2014
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment