How to do registration and resampling of T1W 3D image and the tumor mask with MNI152 template?

Dear forum members,

We have T1W 3D image and its tumor mask of shape 240 x 240 x 155. Now we want to register both image and mask to MNI152 template. MNI152 (mni_icbm152_t1_nlin_asym_09a) shape is 197 x 233 x 189. We are doing linear registration and image is registered but it seems mask is not registered. How do to registration and resampling of the tumor mask so that it will match with the registered image?

We will highly appreciate your help.

Thank you so much.

Thanks & regards,
Sunita

Hello @Sunita_Patel,

If the goal is to transfer the segmentation associated with the T1W image then your code should look something like the following:

fixed_image = mni_template
moving_image = t1w_image
tx = register(fixed_image, moving_image)
resampled_t1w_segmentation = sitk.Resample(t1w_segmentation, fixed_image, tx, sitk.sitkNearestNeighbor)
resampled_t1w = sitk.Resample(t1w_image, fixed_image, tx)

The transformation returned by the registration maps points from the fixed_image coordinate system to the moving_image coordinate system. For additional details see the registration overview documentation.

2 Likes

Hi Ziv,

Thank you so much.
It is working great!

The code we are using is given below. We are using all the default parameters. For the best results which parameters need to be adjusted and anything else we need to include or delete from here?

Could you please help us on how to perform noise removal, and bias correction?

Is the motion correction same as registration?

Thanks again Ziv for your help.

Best regards,
Sunita

    fixed_image =  sitk.ReadImage(img_fixed, sitk.sitkFloat32)
    moving_image = sitk.ReadImage(img_moving, sitk.sitkFloat32)
    moving_mask = sitk.ReadImage(msk_image, sitk.sitkInt8)
    initial_transform = sitk.CenteredTransformInitializer(fixed_image, moving_image, sitk.Euler3DTransform(),
                                                         sitk.CenteredTransformInitializerFilter.GEOMETRY)
    registration_method = sitk.ImageRegistrationMethod()

    # Similarity metric settings.
    registration_method.SetMetricAsMattesMutualInformation(numberOfHistogramBins=50)
    registration_method.SetMetricSamplingStrategy(registration_method.RANDOM)
    registration_method.SetMetricSamplingPercentage(0.40, seed=100)
    registration_method.SetInterpolator(sitk.sitkAffine)
    # Optimizer settings.
    registration_method.SetOptimizerAsGradientDescent(learningRate=0.8, numberOfIterations=60,
                                                      convergenceMinimumValue=1e-6,
                                                      convergenceWindowSize=10)
    registration_method.SetOptimizerScalesFromPhysicalShift()
    # Setup for the multi-resolution framework.
    registration_method.SetShrinkFactorsPerLevel(shrinkFactors=[4, 2, 1])
    registration_method.SetSmoothingSigmasPerLevel(smoothingSigmas=[2, 1, 0])
    registration_method.SmoothingSigmasAreSpecifiedInPhysicalUnitsOn()
    # Don't optimize in-place, we would possibly like to run this cell multiple times.
    registration_method.SetInitialTransform(initial_transform, inPlace=False)

    final_transform = registration_method.Execute(fixed_image, moving_image)
    moving_resampled = sitk.Resample(moving_image, fixed_image, final_transform)
    
    # Resample segmentation
    resampled_t1w_segmentation = sitk.Resample(moving_mask, fixed_image, final_transform, sitk.sitkNearestNeighbor)

Hello @Sunita_Patel,

In terms of parameter settings, you will have to play with those to find the optimal values for your setup (there are no values that are appropriate for all tasks). Also, before evaluating the registration quality, you need to explicitly specify the desired target registration error appropriate for your analysis. The only thing that is certain is that it isn’t an error of zero mm, but it is possible that any value below 2mm is sufficiently good. This bound is task specific.

For bias correction see the N4 bias field correction example.

Motion correction is not registration. Motion correction is a specific task, it can be done using registration. For a general overview see the paper Motion correction in MRI of the brain.

2 Likes

Thank you so much Ziv.

Hi Ziv,

We tried calculating registration error as given in the following code. Is this the correct way of calculating TRE? The way we are taking fixed and moving points are correct or wrong?

Initial transformation we are getting 0 and final transformation we are getting mean(std) as 29(11) mm. Could you please help us why we are getting such as high values. You suggested that it should be less than 2 mm. Right?

Thank you so much for your help Ziv.

Best regards,
Sunita

# (corresponding points in the fixed and moving images).
    fixed_points = ru.generate_random_pointset(image=fixed_image, num_points=100)
    moving_points = [initial_transform.TransformPoint(p) for p in fixed_points]

    (
        initial_errors_mean,
        initial_errors_std,
        initial_errors_min,
        initial_errors_max,
        _,
    ) = ru.registration_errors(
        initial_transform,
        fixed_points,
        moving_points,
        display_errors=False,
    )
    print(
        f"After initialization, errors (TRE) in millimeters, mean(std): {initial_errors_mean:.2f}({initial_errors_std:.2f}), max: {initial_errors_max:.2f}"
    )

    final_errors_mean, final_errors_std, _, final_errors_max, _ = (
        ru.registration_errors(
            final_transform,
            fixed_points,
            moving_points,
            min_err=initial_errors_min,
            max_err=initial_errors_max,
            display_errors=False,
        )
    )


    print(
        f"After registration, errors in millimeters, mean(std): {final_errors_mean:.2f}({final_errors_std:.2f}), max: {final_errors_max:.2f}"
    )

Hello @Sunita_Patel,

“You suggested that it should be less than 2 mm. Right?” - You misunderstood the guidance. You need to identify the required accuracy for your application, 2mm is an example. I do not know what is the requirement for your application it could be 1mm but likely not 10mm. In any case, once that is defined you can judge whether a registration is sufficiently good for your task.

With respect to target registration error and its computation, please see this notebook. Highly recommend that you read the publications referenced there where the fundamentals of registration evaluation are described.

1 Like

I will definitely go through the reference paper and the notebook.
Thanks again Ziv.