Unexpected zero-valued slices when computing Jacobian determinant from SCT warp field

MRI Processing & SCT
1

Hi SCT team,

I am trying to compute a Jacobian determinant image from the deformation field generated by sct_register_to_template.

My registration command was approximately:

sct_register_to_template \
  -i data.nii.gz \
  -s seg.nii.gz \
  -ldisc labels_disc.nii.gz \
  -c t1 \
  -ofolder reg

I then tried to compute the Jacobian determinant using ANTs:

CreateJacobianDeterminantImage \
  3 \
  reg/warp_anat2template.nii.gz \
  spinal_jacobian.nii.gz \
  0 \
  0

The command runs, but the resulting spinal_jacobian.nii.gz contains many zero-valued voxels, including some entire z-slices with value 0. Some of these zero-valued areas are inside the PAM50 template signal range, and the corresponding warp_anat2template.nii.gz voxel has a non-zero value.

For example, in FSLeyes I can see a voxel where:

warp_anat2template: non-zero value
spinal_jacobian: 0.0
PAM50_t1: spinal cord area with non-zero signal

557e39166a9517198e14a67475c49254
557e39166a9517198e14a67475c492541920Ă—1228 149 KB

So I am wondering whether I am using the wrong command or whether the SCT warp field needs to be converted before being passed to CreateJacobianDeterminantImage.

The header of the warp field is:

dim0        5
dim1        141
dim2        141
dim3        991
dim4        1
dim5        3
datatype    16
pixdim0     -1.000000
pixdim1     0.500000
pixdim2     0.500000
pixdim3     0.500000
pixdim4     0.000000
pixdim5     0.000000
intent      Vector
intent_code 1007

PrintHeader also shows:

Spacing [0.5, 0.5, 0.5]
Size : 141 141 991
Image Dimensions : [141, 141, 991]
Intensity Range  : [-100000, 3.40282e+38]
Mean Intensity   : 4.17588e+35

My questions are:

  1. Is CreateJacobianDeterminantImage the correct way to compute the Jacobian determinant from warp_anat2template.nii.gz generated by SCT?
  2. Does the SCT warp field need to be converted from the 5D vector NIfTI format before Jacobian calculation?
  3. Are the extreme values in the warp field, such as -100000 and 3.40282e+38, expected padding values? Could they explain the zero-valued Jacobian slices?
  4. Is there an SCT-recommended way to generate a Jacobian determinant image from sct_register_to_template outputs?

I have attached a minimal set of files for testing via Google Drive:
testing data

Thanks!

Xiaomin

2

Hi Xiaomin,

Thank you so much for your questions. I’m happy to take a look. :blush:

Initial questions

    1. Is CreateJacobianDeterminantImage the correct way to compute the Jacobian determinant from warp_anat2template.nii.gz generated by SCT?
    1. Is there an SCT-recommended way to generate a Jacobian determinant image from sct_register_to_template outputs?

I will have to check in with my colleagues, as I have not used the CreateJacobianDeterminantImage script before.

    1. Does the SCT warp field need to be converted from the 5D vector NIfTI format before Jacobian calculation?

I believe that SCT should output the warping field in the correct 5D format to begin with (see this section of our Warping fields documentation page). I think the 5D format outlined there matches the snippet of the header you shared (i.e. t=1 for dim4):

dim0        5
dim1        141
dim2        141
dim3        991
dim4        1
dim5        3
  1. Are the extreme values in the warp field, such as -100000 and 3.40282e+38, expected padding values? Could they explain the zero-valued Jacobian slices?

Yes, those are expected padding values. I am not sure if they explain the zero-valued slices, but we can put a pin in that for now and come back to it. :slight_smile:

ANTs call

Just for my own understanding and for posterity, I wanted to make sure I understand the positional arguments in the function call you shared using the ANTs documentation. The script usage seems to be:

CreateJacobianDeterminantImage imageDimension deformationField outputImage [doLogJacobian=0] [useGeometric=0] [deformationGradient=0]

And so your command is something like:

CreateJacobianDeterminantImage \
  3 \                              # 1. imageDimension: 3D processing
  reg/warp_anat2template.nii.gz \  # 2. deformationField
  spinal_jacobian.nii.gz \         # 3. outputImage
  0 \                              # 4. doLogJacobian=0
  0 \                              # 5. useGeometric=0 -> standard Jacobian calc, not the geometric method
  0                                # 6. deformationGradient=0 -> implied default value

Aside: I notice that the documentation for ANTs contains the following notice:

Update 2025-05-29 A bug affecting Jacobian determinants in some data is fixed in #1887. The “geometric” Jacobian calculation is not affected. […] If you ran without the geometric option […] your data might be affected. Details in #1884.

Can you double-check whether your version of the function is subject to this bug? Specifically:

A quick way to check without recompiling or upgrading is to re-run CreateJacobianDeterminantImage with the geometric option on. The methods are not equivalent so the results will not be identical, but they should be similar - within the brain, they are mostly within 1% of each other in my tests.

Notes

Assuming that the bug is not the issue, then my first instinct (you have already done this) would be to overlay both the warping field and the output Jacobian image to see if the gaps are also present in the warping field itself (e.g. an issue with the registration process) or if the warping field looks okay but the gaps are introduced by the ANTs script.

My second thought (purely instinctual without a deep understanding of the ANTs script) is that perhaps because the input data is being warped from its native space to a larger space (PAM50), the non-zero slices at specific points correspond to the slices present in the input (with the zero slices corresponding to slices that would otherwise be interpolated as part of the application of the warping field). I don’t know how likely this is given that you describe the warping field as having non-zero values at those slices.

I will try to look more into this when I find a moment. I just wanted to provide an off-the-cuff reply now to let you know I’m looking at this and will follow-up. :slight_smile:

Kind regards,
Joshua

3

Hi Joshua,

Thank you very much for your detailed reply and suggestions.

I checked my ANTs version, and it appears to be the latest version, so I suspect the issue is probably not related to the older ANTs Jacobian determinant bug you mentioned.

I also tried changing the useGeometric option from 0 to 1. However, the result still shows alternating zero-valued slices along the z direction.

To be clear, I am not certain that CreateJacobianDeterminantImage is the correct way to compute the Jacobian determinant from an SCT-generated warp field. My understanding is still limited :baby:, but since SCT registration to the template appears to use ANTs-compatible transformations, and the nonlinear registration step seems to involve BSplineSyN, I thought this ANTs command might be a reasonable thing to try.

For brain MRI data, I have previously used ants.registration.create_jacobian_determinant_image, and the results seemed to work as expected. However, for the spinal cord case, I could not find an official SCT command or documented SCT workflow for extracting the Jacobian determinant from warp_anat2template.nii.gz, so I was experimenting with this approach.

Do you know whether SCT has an official or recommended method for computing the Jacobian determinant from the template-registration warp field? Or is there a preferred way to handle the SCT warp field before computing the Jacobian?

Thanks again for all your help — I really appreciate you taking the time to look into this. :blush:

Best,
Xiaomin

4

Hi @Xiaomin_Lin
This doesn’t answer your question, sorry, (I’m looking into) :upside_down_face:,
but another easy/didactic way to view the warps files (not Jacobians) is through ITKsnap: Multi-component Display/grid

image
image1600Ă—875 278 KB

Best,

Nilser

5

Hi @joshuacwnewton
From the multimodal-registration project, I found this script to calculate jacobians,
I tried to apply it in the @Xiaomin_Lin warps, but I think it only works in multimodal registrations :face_with_diagonal_mouth:, and not to-template registrations.

I’ll keep looking into.

Best,

Nilser

6

Hi @nlaines,

Thank you very much for looking into this and for the helpful suggestion.

I will try visualizing the warp field in ITK-SNAP using the Multi-component Display / grid option to better understand whether the deformation field itself looks reasonable. That sounds like a useful debugging step before focusing on the Jacobian output.

Thank you also for checking the Jacobian script from the multimodal-registration project. I see what you mean that it may be specific to multimodal registrations and may not directly apply to the to-template registration case. If you happen to find a workflow or script that is appropriate for sct_register_to_template outputs, especially warp_anat2template.nii.gz, I would be very grateful.

Thanks again for your help!

Best,
Xiaomin

1 Like
7

Hi @Xiaomin_Lin

I understand that you are interested in computing the Jacobian from sct_register_to_template. This is also something I would like to explore further, but I have not been able to generate it directly in a fully consistent way so far.

I have been thinking about it, and here is what I have found:

  1. @joshuacwnewton when I mentioned that this script only worked with sct_register_multimodal registrations, I was referring to previous experiments in multimodal settings where I obtained reasonable-looking Jacobian maps. However, when applying the same approach to sct_register_to_template warps, it did not behave as expected, and the same issue appeared when using CreateJacobianDeterminantImage (as @Xiaomin_Lin mentions)

  2. I explored an alternative approach that may or may not be strictly correct, but it appears to work in practice. Since the Jacobian reflects deformation introduced by nonlinear registration, it is generally more meaningful in nonlinear steps such as BSplineSyN. Based on this idea, I looked at splitting the standard sct_register_to_template pipeline.

The standard usage of sct_register_to_template is:

sct_register_to_template \
  -i t1_spinal.nii.gz \
  -s t1_spinal_seg.nii.gz \
  -ldisc t1_spinal_seg_labeled_discs.nii.gz \
  -ofolder reg

Internally, this performs two main steps:

  • step=1,type=imseg,algo=centermassrot,metric=MeanSquares:\
  • step=2,type=seg,algo=bsplinesyn,metric=MeanSquares,iter=3

My idea is to explicitly separate these two steps by adjusting -param.

First step (initial alignment):

sct_register_to_template \
  -i t1_spinal.nii.gz \
  -s t1_spinal_seg.nii.gz \
  -ldisc t1_spinal_seg_labeled_discs.nii.gz \
  -param \
step=1,type=imseg,algo=centermassrot,metric=MeanSquares \
  -ofolder pam2steps

Before running the second step, we need to bring the SC segmentation into
the intermediate space obtained after step 1. For that, we apply the
transformation:

sct_apply_transfo \
  -i t1_spinal_seg.nii.gz \
  -d /Path/PAM50/template/PAM50_t1.nii.gz \
  -w pam2steps/warp_anat2template.nii.gz \
  -o pam2steps/anat2template_seg_from-applytransfo.nii.gz \
  -x nn

Then we run the second step (using sct_register_multimodal):

sct_register_multimodal \
  -i pam2steps/anat2template.nii.gz \
  -d /Path/PAM50/template/PAM50_t1.nii.gz \
  -iseg pam2steps/anat2template_seg_from-applytransfo.nii.gz \
  -dseg /Path/PAM50/template/PAM50_cord.nii.gz \
  -param \
step=1,type=seg,algo=bsplinesyn,metric=MeanSquares,iter=3 \
  -o pam2steps/t1_second_stage_2_pam.nii.gz

Finally, we can compute the Jacobian only from the second step:

CreateJacobianDeterminantImage 3 \
pam2steps/warp_anat2template2PAM50_t1.nii.gz \
pam2steps/jacobian_t1_second_stage_to_pam.nii.gz \
1 0

So in this approach, we are effectively isolating the Jacobian
corresponding only to the BSplineSyN deformation.

image
image1853Ă—961 148 KB

I am not completely sure whether this is the correct way to interpret it
for the effect we want to observe, but I hope it is helpful :slight_smile:

Best,

Nilser

8

Thank you so much for your input, @nlaines!

Since the Jacobian reflects deformation introduced by nonlinear registration, it is generally more meaningful in nonlinear steps such as BSplineSyN

In that case, I am wondering if it is not so much “sct_register_to_template vs. sct_register_multimodal” then, but “nonlinear vs. linear registration” steps?

In other words, the Jacobian image is only valid for registrations using particular -param settings, and perhaps the default -param settings differ between the two scripts (or that perhaps different -param settings were used in different test runs) causing the unexpected output.

I think there may be a way to recover the intermediate warping fields without splitting up the command? I would try to run -r 0 (turning off the removal of temporary files) then looking in /tmp for the most recent SCT-generated tmpdir, then finding the warping fields for each step. (I think that SCT builds warping fields for each step individually, then concatenates them together at the very end.)

I am wondering if this is a big enough motivation to expose the intermediate warping fields in the output as well rather than only the concatenated output. That could make it easier to understand the impact of each individual step, and aid in debugging in the future. :thinking:

Certainly something to test. I can bring it up at the next SCT developer meeting, too.

1 Like
9

Hi @nlaines and @joshuacwnewton ,

Thank you both so much for the detailed explanations and suggestions! This was really helpful, and I feel like I now have a much better understanding of what might be going on, especially regarding the difference between the full concatenated warp field and the intermediate nonlinear registration step.

Nilser, your workflow makes a lot of sense to me. I’ll try separating the initial alignment from the BSplineSyN step and compute the Jacobian only from the nonlinear warp field.

Joshua, thanks also for your help! It would be great to see SCT support this kind of Jacobian-related analysis more directly in the future. I have come across some studies suggesting that Jacobian determinants may contain useful information for assessing atrophy in patient groups, so I think having a clearer and more standardized way to compute and interpret them in SCT would be very valuable.

Really appreciate your help!

Best,
Xiaomin

2 Likes
10

Thank you for the feature request! I will note this down in an issue in the SCT repository, and bring it up in our next SCT developer meeting on Tuesday. It could certainly be an interesting project for us to work on. :slight_smile:

I appreciate your involvement on the forum and your interest in improving SCT! :hearts:

Kind regards,
Joshua

11

For the sake of adding additional context to the newly-created issue (Add built-in support for Jacobian determinant image analysis on SCT warping fields · Issue #5272 · spinalcordtoolbox/spinalcordtoolbox · GitHub), would you happen to have any links you came across during your literature review on hand as a jumping off point for further research? :slight_smile:

12

Hi Joshua @joshuacwnewton ,

Sorry for the delayed reply—I completely missed your message.

Yes, during my literature review I came across several papers where Jacobian determinants (or Jacobian maps derived from nonlinear registrations) are used as quantitative measures of local tissue expansion/atrophy - most of which focus on the brain. A few references that might be useful as a starting point are:

  1. Leow et al., 2007, Statistical properties of Jacobian maps and the realization of unbiased large-deformation nonlinear image registration (IEEE Transactions on Medical Imaging, 26(6):822–832).
    • This is probably one of the foundational paper. It discusses the statistical properties of Jacobian maps and why unbiased nonlinear registration is important when using Jacobian determinants as morphometric measures.
  2. Xu et al., Nature (2023) – Normative spatiotemporal fetal brain maturation with satisfactory development at 2 years.
    Normative spatiotemporal fetal brain maturation with satisfactory development at 2 years | Nature
    • They use tensor-based morphometry (TBM), where Jacobian determinants derived from nonlinear registration are analyzed to characterize local brain growth and developmental changes.
  3. Brain (2020) – Sex effects on brain structure in de novo Parkinson’s disease: a multimodal neuroimaging study.
    https://academic.oup.com/brain/article/143/10/3052/5912110
    • This study uses deformation-based morphometry (DBM) to quantify regional atrophy patterns in Parkinson’s disease, with Jacobian-derived measures serving as the primary morphometric metric.
  4. Elife (2015) – Network structure of brain atrophy in de novo Parkinson’s disease
    Network structure of brain atrophy in de novo Parkinson's disease | eLife
    • Another interesting study used DBM
  5. “Mapping pediatric spinal cord development with age” (2022).
    Mapping Pediatric Spinal Cord Development with Age - PMC
    • This one may be particularly relevant to the SCT issue because it applies Jacobian-based morphometry directly to the spinal cord. They computed Jacobian determinant maps from ANTs registrations to characterize normal spinal cord growth across vertebral levels. While the work focuses on normative development rather than disease-related atrophy, it demonstrates that Jacobian-based morphometry is feasible in the spinal cord and can capture spatially varying tissue changes.

Building on that, I think it would be interesting to investigate whether Jacobian determinant maps could also provide useful information in pathological spinal cord datasets, for example as a complementary voxel-wise measure of local tissue expansion or atrophy. Although SCT already provides established morphometric metrics such as CSA, a Jacobian-based workflow could potentially reveal spatial patterns of change that are not captured by regional summary measures alone :thinking:.

Best,
Xiaomin

1 Like
13

Those extreme values you were asking about will probably look different once you pull the intermediate BSplineSyN warp—concatenated deformation fields tend to have weird stuff like that from combining multiple registration stages.