ANTsX / ANTsR

R interface to the ANTs biomedical image processing library
https://antsx.github.io/ANTsR
Apache License 2.0
127 stars 35 forks source link

What is #352

Closed muratmaga closed 2 years ago

muratmaga commented 2 years ago

I am testing a two step registration, where I do a SyN registration with lower resolution dataset and I use the output from that to initialize the same registration with the original images.

syn = antsRegistration(fixed=tmplt, 
                             moving=mov.img, 
                             typeofTransform = "SyN", 
                             initialTransform = lm.tx$transform)

#same registration at higher res
      syn2 = antsRegistration(fixed=tmplt2, 
                              moving=mov.img2, 
                              typeofTransform = "SyN", 
                              initialTransform = syn$fwdtransforms) 

Frequently, but not for all samples, I am getting this message at the end of syn2 step:

transform composite list may not be invertible - return all transforms and leave it to user to figure it out

I am not entirely sure what I am supposed to do with this information. Does this mean registration hasn't worked as intended? Should I ignore the result of these? Or this is harmless?

If it matters, my final step is to apply these to a point set on the template image to get them on the moving image:

syn2LM = antsApplyTransformsToPoints(dim = 3, 
                                          points = ref.lm, 
                                          transformlist = syn2$fwdtransforms)
ntustison commented 2 years ago

As long as you're keeping track of the transforms and understand the issues involved, you can ignore the message. The issue is that in supplying the SyN stage 1 forward transform as the initial transform of your second stage, the registration program is not simultaneously supplied with the corresponding inverse transform (and is not even given the information that the supplied transform is invertible). So that warning message is given just to give notice to the user.

muratmaga commented 2 years ago

Thanks Nick. How can I pass the inverse transform and/or specify that the provided transform is invertable? I will need the accurate inverse transform from second syn at some point.

ntustison commented 2 years ago

You don't. It's not necessary for the registration. As I mentioned above, the user is responsible for keeping track of the transforms in these types of situations. After the second stage, you can use antsApplyTransforms to concatenate transforms, if necessary.

muratmaga commented 2 years ago

One more clarification request.

Is what you describe (that I need to keep track and provide the correct sequence of transforms if I need to use two stage registration like this) specific only to inverse transform or both fwd and inverse transform of second syn (asking strictly for my example above)? Because I see the fwdtransforms seem to have the correct number (and the order) of transforms as far as I can tell), but not the inverse.

$fwdtransforms
[1] "/scratch/RtmpHewqup/file833e4a547b863Warp.nii.gz"
[2] "/scratch/RtmpHewqup/file833e4a547b862GenericAffine.mat"
[3] "/scratch/RtmpHewqup/file833e4a547b861Warp.nii.gz"
[4] "/scratch/RtmpHewqup/file833e4a547b860GenericAffine.mat"

$invtransforms
[1] "/scratch/RtmpHewqup/file833e4a547b860GenericAffine.mat"
[2] "/scratch/RtmpHewqup/file833e4a547b861Warp.nii.gz"
[3] "/scratch/RtmpHewqup/file833e4a547b862GenericAffine.mat"
[4] "/scratch/RtmpHewqup/file833e4a547b863InverseWarp.nii.gz"
[5] "/scratch/RtmpHewqup/file833e4a547b863Warp.nii.gz" 

Thanks again...

ntustison commented 2 years ago

Yes, that's exactly right. Because you specified the initial forward transform in the antsRegistration call, that actually gets added to the returned transform list.

muratmaga commented 2 years ago

Awesome, thanks for the great explanation.

muratmaga commented 2 years ago

I am sorry to reopen this, but I now need the inverse transform (fwd transform works fine), but I cannot seem to figure out the correct order of the transform to make it work. I created a small example using built in data. I appreciate if you can help me find out the correct combination:


fi <- antsImageRead(getANTsRData("r16") )
mi <- antsImageRead(getANTsRData("r64") )
fmask = getMask(fi)
mmask = getMask(mi)

fi2 = resampleImage(fi, c(2,2))
mi2 = resampleImage(mi, c(2,2))

syn1 = antsRegistration(fi2, mi2, "SyN")
syn2 = antsRegistration(fi, mi, "SyNCC", syn1$fwdtransforms)

#this works. Just a sanity check
test.mask = antsApplyTransforms(fixed=fi, moving=mmask, transformlist = syn2$fwdtransforms, 'genericLabel')
plot(fi, test.mask, alpha=.8)

#cant figure out the correct order of transforms for inverse
test.mask = antsApplyTransforms(fixed=mi, moving=fmask, transformlist = c(syn1$invtransforms, syn2$invtransforms[3:4]), 'genericLabel')
plot(mi, test.mask, alpha=.8)
ntustison commented 2 years ago

First thing to note is that your displacement fields for syn1 are sampled to the space of the fixed image for that registration pair so you're going to need to resample them to your second fixed image if you want to carry them through with the inclusion of syn2 transforms. So taking your example:

fi <- antsImageRead(getANTsRData("r16") )
mi <- antsImageRead(getANTsRData("r64") )
fmask = getMask(fi)
mmask = getMask(mi)

fi2 = resampleImage(fi, c(2,2))
mi2 = resampleImage(mi, c(2,2))

syn1 = antsRegistration(fi2, mi2, "SyN")
syn2 = antsRegistration(fi, mi, "SyNCC", syn1$fwdtransforms)

# Resampling --- we'll do both directions for completeness

# (Note file names simply used for convenience here.  Need to change locally.)
inv_warp_file <- "/Users/ntustison/syn1_invwarp_resampled.nii.gz"
warp_file <- "/Users/ntustison/syn1_warp_resampled.nii.gz"

syn1_invwarp = antsImageRead( syn1$invtransforms[2] )
syn1_invwarp_split <- splitChannels( syn1_invwarp )
syn1_invwarp_resampled <- list()
syn1_invwarp_resampled[[1]] <- resampleImageToTarget( syn1_invwarp_split[[1]], fi )
syn1_invwarp_resampled[[2]] <- resampleImageToTarget( syn1_invwarp_split[[2]], fi )
syn1_invwarp_resampled <- mergeChannels( imageList = syn1_invwarp_resampled )
antsImageWrite( syn1_invwarp_resampled, inv_warp_file )

syn1_warp = antsImageRead( syn1$fwdtransforms[1] )
syn1_warp_split <- splitChannels( syn1_warp )
syn1_warp_resampled <- list()
syn1_warp_resampled[[1]] <- resampleImageToTarget( syn1_warp_split[[1]], fi )
syn1_warp_resampled[[2]] <- resampleImageToTarget( syn1_warp_split[[2]], fi )
syn1_warp_resampled <- mergeChannels( imageList = syn1_warp_resampled )
antsImageWrite( syn1_warp_resampled, warp_file )

# Create a combined set of transforms forward and backward

synAll <- list()
synAll$fwdtransforms <- c( syn2$fwdtransforms[1], syn2$fwdtransforms[2], warp_file, syn1$fwdtransforms[2] )
synAll$invtransforms <- c( syn1$invtransforms[1], inv_warp_file, syn2$invtransforms[3], syn2$invtransforms[4] )

# Test

# This checks out locally
test.mask = antsApplyTransforms( fixed = fi, moving = mmask, transformlist = synAll$fwdtransforms )
plot( fi, test.mask, alpha = 0.8 )

# This checks out locally
#  IMPORTANT:  Don't forget the whichtoinvert option
test.mask = antsApplyTransforms( fixed = mi, moving = fmask, transformlist = synAll$invtransforms, whichtoinvert = c( TRUE, FALSE, TRUE, FALSE ) )
plot( mi, test.mask, alpha = 0.8 )
muratmaga commented 2 years ago

Awesome, thank you!