Parameter Maps

In the previous section we saw how to configure an entire multi-resolution registration process with a single parameter. Here, we will examine how SimpleElastix configures registration components internally. Knowing the internal mechanisms is not strictly necessary to register simple anatomical structures, but it will help you solve complex problems that require problem-specific tuning.

Elastix introduces the concept of a parameter map to configure the registration procedure. A parameter map is a collection of key-value pairs that atomically defines the components of the registration and any settings they might require. Only input images and output options need to be specified seperately. Elastix will read a given parameter map and load the specified components at runtime.

The original elastix and transformix command line programs read text files from disk in which parameters are specified according to the following format:

(ParameterName ParameterValue)

The parameter may be either single- or vector-valued. If the value is of type string then the value will be quoted:

(ParameterName "value1" ... "valueN")

If the value is of numerical type it will be unquoted:

(ParameterName 123)

SimpleElastix can read these types of files as well but further introduces native data objects for parameter files in all target languages. This means that instead of editing text files you can programmatically configure registration components. For example, we could initialize a parameter object and start populating it with values:

import SimpleITK as sitk

p = sitk.ParameterMap()
p['Registration'] = ['MultiResolutionRegistration']
p['Transform'] = ['TranslationTransform']
...

and so on. ElastixImageFilter has default settings that allows us to get started right away.

The Default Parameter Maps

In the Hello World example we obtained a registered image by running

resultImage = sitk.Elastix(sitk.ReadImage('fixedImage.nii'),  \
                           sitk.ReadImage('movingImage.nii'), \
                           'translation')

Internally, ElastixImageFilter uses a pre-configured parameter map to register images. The parameter map is obtained by an internal call to sitk.GetDefaultParameterFile('translation'). This function provides parameter maps for rigid, affine, non-rigid and groupwise registration methods (in order of increasing complexity).

Tip

ElastixImageFilter will register our images with a translation -> affine -> b-spline multi-resolution approach by default. We simply leave out the call to SetParameterMap to achieve this functionality.

import SimpleITK as sitk

# Functional interface
resultImage = sitk.Elastix(sitk.ReadImage('fixedImage.nii'), sitk.ReadImage('movingImage.nii'))

# Object oritented interface
elastixImageFilter = sitk.ElastixImageFilter()
elastixImageFilter.SetFixedImage(sitk.ReadImage('fixedImage.nii'))
elastixImageFilter.SetMovingImage(sitk.ReadImage('movingImage.nii'))
resultImage = elastixImageFilter.Execute()

We can also retrieve the parameter map ourselves and reconfigure it before passing it back to ElastixImageFilter, allowing us to quickly optimize a registration method for a particular problem:

parameterMap = sitk.GetDefaultParameterMap('translation')

# Use a non-rigid transform instead of a translation transform
parameterMap['Transform'] = ['BSplineTransform']

# Because of the increased complexity of the b-spline transform,
# it is a good idea to run the registration a little longer to
# ensure convergence
parameterMap['MaximumNumberOfIterations'] = ['512']

resultImage = sitk.Elastix(sitk.ReadImage('fixedImage.nii'),  \
                           sitk.ReadImage('movingImage.nii'), \
                           parameterMap)

Tip

We can print parameter maps to console like this:

import SimpleITK as sitk
elastixImageFilter = sitk.ElastixImageFilter()
elastixImageFilter.PrintParameterMap()

We will study other parameter maps more closely in later examples. For now, we simply print the translation parameter map to console and examine its contents.

>>> sitk.PrintParameterMap(sitk.GetDefaultParameterMap("translation"))
ParameterMap 0:
  (AutomaticParameterEstimation "true")
  (CheckNumberOfSamples "true")
  (DefaultPixelValue 0)
  (FinalBSplineInterpolationOrder 2)
  (FixedImagePyramid "FixedSmoothingImagePyramid")
  (FixedImagePyramidSchedule 8 8 8 4 4 4 2 2 2 1 1 1)
  (ImageSampler "RandomCoordinate")
  (Interpolator "LinearInterpolator")
  (MaximumNumberOfIterations 32)
  (MaximumNumberOfSamplingAttempts 8)
  (Metric "AdvancedMattesMutualInformation")
  (MovingImagePyramid "MovingSmoothingImagePyramid")
  (MovingImagePyramidSchedule 8 8 8 4 4 4 2 2 2 1 1 1)
  (NewSamplesEveryIteration "true")
  (NumberOfResolutions 4)
  (NumberOfSamplesForExactGradient 4096)
  (NumberOfSpatialSamples 4096)
  (Optimizer "AdaptiveStochasticGradientDescent")
  (Registration "MultiResolutionRegistration")
  (ResampleInterpolator "FinalBSplineInterpolator")
  (Resampler "DefaultResampler")
  (Transform "TranslationTransform")
  (WriteResultImage "true")

The first thing to note is that the parameter map is enumerated. ElastixImageFilter can take a vector of parameter maps and apply the corresponding registrations sequentially. The resulting transform is called a composite transform since the final transformation is a composition of sequentially applied deformation fields. For example, a non-rigid registration is often initialized with an affine transformation (translation, scale, rotation, shearing) to bring the objects into rough alignment. This makes the registration less suscetible to local minima. We can also ask SimpleElastix to add the individual deformation fields and apply them in one go (but make sure you know what you are doing before opting for this apprach).

Tip

We can add mulitple parameter maps to SimpleElastix like this:

import SimpleITK as sitk
elastixImageFilter = sitk.ElastixImageFilter()
elastixImageFilter.SetParameterMap(sitk.GetDefaultParameterMap('translation'))
elastixImageFilter.AddParameterMap(sitk.GetDefaultParameterMap('affine'))

Note that the first call is a Set method. This deletes any prevously set parameter maps. Subsequent calls to AddParameterMap appends parameter maps to the internal list of parameter maps

Let’s examine the parameters above in detail.

Important Parameters

Registration is the top-level parameter which in this case has been set to MultiResolutionRegistration. A multi-resolution pyramid strategy improves the capture range and robustness of the registration. We will almost always want to use multiple resolutions unless your problem is particularly simple. The basic idea is to first estimate T(x) on a low resolution version of the images and then propagate the estimated deformation to higher resolutions. This makes the registration initially focus on larger structures (the skull and brain hemispheres etc), before focusing on high-frequency information (brain subregions etc) which contain more local minima. FixedImagePyramid, FixedImagePyramidSchedule, MovingImagePyramid, MovingImagePyramidSchedule, and NumberOfResolutions controls the pyramid strategy.

The Transform parameter is set to TranslationTransform which it is optimized with an AdaptiveStochasticGradientDescent optimizer (Klein et al. 2009). SimpleElastix will use this optimizer together with the AdvancedMattesMutualInformation metric by default since this combination work well for a broad range of problems whether mono-modal or multi-modal.

Image intensities are sampled using an ImageSampler, Interpolator and ResampleInterpolator. The sampler is responsible for selecting points in the image to sample. The RandomCoordinate simply selects random positions. The interpolator is responsible for interpolating off-grid posititions during optimization. The LinearInterpolator used here is very fast and uses very little memory.

A BSplineInterpolator of order 2 is used to resample the result image from the moving image once the final transformation has been found. This is a one-time step so the additional computational complexity is worth the trade-off for higher image quality.

Another important parameter is AutomaticParameterEstimation which controls whether the AdaptiveStochasticGradientDescent optimizer should estimate its own convergence parameters or allow you to set them. Automatically obtained parameters work well in most cases and facilitates a complete hands-off approach which is highly recommended. Optimizers can be tricky to tune by hand.

DefaultPixelValue sets value of pixels outside the moving image grid. The rest of the key-value pairs are component specific parameters. There are multiple choices available for each type of component. For example, you can construct an image pyramid with recursive sampling or via Gaussian Smoothing. Each choice has its own pros and cons. Consult the Registration Components secton for a description of all types of available components.