What do itk-snap registrations compute exactly?

[alexisthual]

I am trying to register an T1 anatomy IMG2 acquired at the beginning of an fMRI session with a high-resolution T1 anatomy IMG1 that I used with freesurfer’s reconall to extract surfaces.

I went through the following steps:

1. registered IMG2 on IMG1 in itk-snap (see image B)
2. exported computed registration in itk format
3. transformed itk registration using c3d_affine_tool (see command C below), which yielded the following file:

〉 cat registration.mat
0.965521 0.0573185 0.0413778 29.7438
0.0347123 -0.0362805 -0.867828 176.751
-0.0535802 0.973744 -0.00314146 23.563
0 0 0 1

Finally, I applied the computed registration using the following script in Python:

import nibabel as nib
import numpy as np

anat_high_img = nib.load("img1.nii.gz")
anat_sess_img = nib.load("img2.nii.gz")
registration = np.loadtxt("registration.mat", delimiter=" ")

full_registration = anat_high_img.affine @ registration

registered_anat_sess_img = nib.nifti1.Nifti1Image(
    anat_sess_img.get_fdata(),
    full_registration
)
registered_anat_sess_img.to_filename("img2_registered.nii")

but the image saved is not coregistered to the initial T1 image IMG1 (see picture D below).

I am surprised about this because I used this procedure previously and it had given me co-registered images. I tried a lot of combinations of the registration matrix and the affine functions (using that of anat_sess_img as well, in formulas like full_registration = anat_high_img.affine @ registration @ anat_sess_img.affine), but nothing yielded co-registered images.

Am I missing something here?
What is ITK-snap computing exactly?

Thanks a lot in advance for your help!

A. IMG1 (grayscale colormap) and IMG2 (heat colormap) in freeview

B. Registring IMG1 on IMG2 in itk-snap

C. c3d_affine_tool command

c3d_affine_tool -ref img1.nii.gz -src img2.nii.gz -itk registration.txt -ras2fsl -o registration.mat     

D. Co-registered anatomies are mis-aligned

[jsein]

Hi, I used to do the same procedure, i.e. compute the transformation in ITKSNAP and apply it in FSL.

The main difference I see is that I was exporting the ITKSNAP transformation in C3D format: “Convert3D Transform Files” . That may be the origin of your problems here.

The c3d_affine_tool command is used to convert the C3D transform format to FSL format.

Could you try the exact same thing but choosing “Convert3D Transform Files” instead of “ITK Transform Files” when saving the transformation file in ITKSNAP?

[alexisthual]

Thank you for your answer @jsein !
I had already tried using both formats actually, but it didn’t change anything.

For the record, here are the commands I used, and the input / output files:

ITK Transform File

c3d_affine_tool -ref img1.nii.gz -src img2.nii.gz -itk registration_itk.txt -ras2fsl -o registration_from_itk.mat     

Input

〉 cat registration_itk.txt 
#Insight Transform File V1.0
#Transform 0
Transform: MatrixOffsetTransformBase_double_3_3
Parameters: 1.0247199238239288 0.05537982422696509 -0.04231539935602776 -0.05592859613724214 1.0226895032184466 -0.0014128956388533966 0.04368502987618375 0.04532678510345843 1.151493955202453 4.549996774855572 28.15756508420334 -14.1719031661072
FixedParameters: 0 0 0

Output

〉 cat registration_from_itk.mat 
0.9714 -0.0356308 0.0541816 -100.712
-0.0530699 0.00314632 0.974801 -26.6025
-0.0389417 -0.866961 0.036316 145.833
0 0 0 1

C3D Transform file

c3d_affine_tool -ref img1.nii.gz -src img2.nii.gz registration_c3d.mat -ras2fsl -o registration_from_c3d.mat     

Input

〉 cat registration_c3d.mat
1.02472 0.0553798 0.0423154 -4.55
-0.0559286 1.02269 0.0014129 -28.1576
-0.043685 -0.0453268 1.15149 -14.1719
0 0 0 1

Output

〉 cat registration_from_c3d.mat 
0.9714 -0.0356309 0.0541816 -100.712
-0.0530699 0.00314633 0.9748 -26.6024
-0.0389418 -0.866964 0.0363161 145.834
0 0 0 1

[jsein]

Indeed your two commands with c3d_affine_tool produce a very similar output.
What I do then is to use FSL to apply these transformation matrices produced by c3d_affine_tool:

For your that would give:

flirt -in img2.nii.gz -ref img1.nii.gz -applyxfm -init registration_from_c3d.mat -out img2_2_img1_FSL_spline -interp spline

Could you try this command and see if the same misregistration occurs?

[alexisthual]

It sounds like a very sensible thing to do, but it seems this is still not the right way to go
What is surprising though is that the generated images are different from the ones I showed previously.

Here is a screenshot with img1 in grayscale, and 2 registered images obtained with flirt: the one with the heat colormap was generated using registration_from_c3d.mat, the one with the jet colormap was generated using registration_c3d.mat

[jsein]

I am surprised. I did some tests on a set of images I got and it works on my side.
Coud you double-check that you didn’t modify the affine transform in the header of your moving image before applying the calculated transformation to it?

If you don’t mind sharing a set of your images as a toy example, I could try on my side to see if I see the same behavior.

[alexisthual]

Hummm, you might be right, I just noticed that the registration matrices in my second message are different from the one I reported in the first message. I’ll do another pass.
In the meantime, here is a wetransfer link to IMG1 and IMG2: https://we.tl/t-TqcZlkacXJ

[alexisthual]

You were right @jsein , I got one of the previous steps wrong and using flirt built the correct image:

flirt -in img2.nii.gz -ref img1.nii.gz -applyxfm -init img2_to_img1_itk_fsl.mat -out img2_to_img1_fsl_spline -interp spline

Before I close this topic, I would still like to understand how to reproduce this in Python (because only then will I have understood what is really going on haha).
Could it be that flirt implicitly uses an affine matrix which I’m not aware of?

[jsein]

This website should answer your questions:

Visualization of affine resamplings

[alexisthual]

Hummm, I think something is still missing
I couldn’t “manually” find a combination of products of affine matrices (the registration computed with itk snap and img1 & img2’ affine matrices) that would generate the correct transformation, so I tried all of them.

For each combination of product of affine matrices, I generated a Nifti image, and loaded it in freeview and fsleyes to check if it was aligned with img1. In short: none of them were.
I also tested inverting all or some of these affine matrices when computing their product.

Could it be that flirt uses data I am not aware of? Are there pieces of information in the headers of these images that encode for some translation which would not be present in the affine matrix?

Code snippet to generate all combinations of products:

affines = [registration, img1.affine, img2.affine]

affine_combinations = [
    list(itertools.permutations(np.arange(len(affines)), 1)),
    list(itertools.permutations(np.arange(len(affines)), 2)),
    list(itertools.permutations(np.arange(len(affines)), 3)),
]

inverse_combinations = [
    list(itertools.product([0, 1], repeat=1)),
    list(itertools.product([0, 1], repeat=2)),
    list(itertools.product([0, 1], repeat=3)),
]

affine_inverse_combinations = list(
    itertools.chain(
        *[
            list(
                itertools.product(
                    affine_combinations[i], inverse_combinations[i]
                )
            )
            for i in range(len(affine_combinations))
        ]
    )
)


def get_affine(affine_indices, should_inverse):
    acc = np.eye(4)
    for i in range(len(affine_indices)):
        if should_inverse[i] == 1:
            acc = np.linalg.inv(affines[affine_indices[i]]) @ acc
        else:
            acc = affines[affine_indices[i]] @ acc

    return acc


all_affines = np.array(
    [
        get_affine(affines, should_inverse)
        for affines, should_inverse in tqdm(affine_inverse_combinations)
    ]
)

all_images = [
    nib.nifti1.Nifti1Image(
        img2.get_fdata(), affine
    )
    for affine in all_affines
]

for i, img in enumerate(all_images):
    img.to_filename(
        "/home/alexis/singbrain/data/tmp/registered_imgs/"
        f"img2_to_img1_{affine_inverse_combinations[i]}.nii.gz"
    )

Examples of 78 affine matrices generated by this script:

• img1.affine @ registration
• registration @ img1.affine
• np.linalg.inv(img1.affine) @ registration @ img2.affine

[jsein]

I am not familiar with the code you expose here. It is always a bit tricky to move from one software to another, as for instance some software apply the transformation to the sform matrice while other touch only the qform matrice, some copy sform to the form or vice-versa.

If you stay with FSL, and FSLeyes for visualization, you can play with the “nudge” tool where you can manually set some affine transformation, see how the image moves and see either the result on the image affine or save the transformation matrice. Maybe that will help to understand the transformation process?

[alexisthual]

Sooo, I think I finally understood what is going on

Say in ITK-Snap you:

• first load an IMG1, then an IMG2
• compute a registration of IMG2 onto IMG1
• export it in C3D format as R.mat

This actually yields a matrix R which describes an affine transformation to be applied to go from the millimeter space of IMG1 to that of IMG2, as described in this figure:

Therefore, here is a simple Python snippet to co-register IMG2 with IMG1:

import nibabel as nib
import numpy as np

R = np.loadtxt("R.mat", delimiter=" ")
img2 = nib.load("img2.nii.gz")

img2_registered_on_img1 = nib.nifti1.Nifti1Image(
    img2.get_fdata(),
    np.linalg.inv(R) @ img2.affine
)

img2_registered_on_img1.to_filename(
    "img2_registered_on_img1.nii.gz"
)

Loading IMG1 and the registered IMG2 in Freesurfer or FSLeyes allows to check that they are aligned.

Here are a few takeaways:

• contrary to what the interface suggests in ITK Snap (in my case, in the registration panel, it says “Moving registration layer: IMG2”), the exported registration actually describes a transformation from IMG1 to IMG2 ; one should invert this matrix if they want to co-register IMG2 onto IMG1
• using c3d_affine_tool -ras2fsl does something tricky: not only does it change the origin of the millimeter space to fit that of FSL, it also inverts the input affine transformation
• this is why flirt works out of the box with affine matrices outputted by c3d_affine_tool. However, I would advocate this should be considered non-standard. In particular, if you use flirt with registration matrices that were not parsed with c3d_affine_tool, the outputted transformation of IMG2 won’t be aligned to the reference image IMG1
• note that the described procedure does not resample IMG2 onto IMG1 (ie their number of voxels can be different), contrary to what flirt does

My humble opinion is that a “generic” pipeline should not rely on the FSL space, and therefore not use either c3d_affine_tool or flirt. To conclude, in the future, I’ll use the the itksnap/python procedure described in this post when I need to manually co-register images!

Thank you a lot for your time and help @jsein !

[jsein]

Thank you @alexisthual for this clear explanation! Good to know!

[TomW]

Hi @alexisthual,

thanks a lot, this thread really helped me!

I am using xmf.txt files (ITK, fMRIprep output) to bring additional files in the tw1 space.
For that I am

1. converting the .txt files with c3d_affine_tool to a mat file and then
2. apply it with a similar python snippet as in your post.

Initially I hoped to apply the xmf.txt files without conversion but it didn’t work.
I was wondering what happens during the conversion and was hoping for your insight:

The xmf file looks like this:

#Insight Transform File V1.0
#Transform 0
Transform: AffineTransform_float_3_3
Parameters: 0.999843 0.00639795 0.0164898 -0.00433495 0.992491 -0.122236 -0.017148 0.122145 0.992364 1.76259 0.532735 -1.19536
FixedParameters: 0 0 0

Then, after the conversion with the c3d_affine_tool, the .mat file looks like this:

command:

c3d_affine_tool -itk img_ses-2_from-orig_to-T1w_mode-image_xfm.txt \
    -o ses-2_to_t1w.mat 

output:

[[ 0.999843    0.00639795 -0.0164898  -1.76259   ]
 [-0.00433495  0.992491    0.122236   -0.532735  ]
 [ 0.017148   -0.122145    0.992364   -1.19536   ]
 [ 0.          0.          0.          1.        ]]

c3d_affine_tool seems to have flipped the signs of the x and y axes of the rotation matrix.
Is that because of an implicit change from LAS to RAS convention?

As you see I am not offering the tool any information besides the xmf file, and also no -ras2fsl flag.
Why does c3d_affine_tool sill convert something by default?

Thanks in advance!

Cheers
Tom

[alexisthual]

Hi @TomW !

I’m a bit hazy on the details now, but I think that if you want to use my python snippet above then you should probably not use c3d_affine_tool (or flirt), and instead work directly with the output of ITK-Snap.

c3d_affine_tool seems to have flipped the signs of the x and y axes of the rotation matrix.
Is that because of an implicit change from LAS to RAS convention?

Indeed, I make a somewhat similar point in the 2nd bullet point of my “takeaways” above, although I don’t really understand what is going on. As a rule of thumb I would try to use as few of these tools as possible because it feels to me that they very often implement different defaults that cannot be easily setup together by a beginner like me!

Good luck

[TomW]

Hey @alexisthual

Thanks for the reply!
The registration looks fine now so I believe it works well, even though I don’t have a clear understanding of all underlying processes.

Yes, sticking to one software seems reasonable.
However, if I read the documentation right, fmriprep computes the transformations from individual sessions to the t1w-space (in case of multiple anatomical images) with freesurfers “robust_template”.
But the output is in ITK format, which already confuses me in terms of the right software to use.

Anyways, thanks again and best wishes