Image Registration Using the 'NiftyReg' Library

Provides an 'R' interface to the 'NiftyReg' image registration tools < http://sourceforge.net/projects/niftyreg/>. Linear and nonlinear registration are supported, in two and three dimensions.


Build Status Build status

The RNiftyReg package is an R-native interface to the NiftyReg image registration library developed within the Translational Imaging Group at University College London. The package incorporates the library, so it does not need to be installed separately, and it replaces the NiftyReg command-line front-end with a direct, in-memory bridge to R, based on Rcpp.

This README file primarily covers version 2.0.0 of the package and later. The interface was substantially reworked in that version to make it more natural and less verbose, and earlier versions are incompatible. Information on moving from prior versions of RNiftyReg to 2.x is included at the end of this file.

The package can be installed from CRAN, or the latest development version obtained directly from GitHub:

## install.packages("devtools")
devtools::install_github("jonclayden/RNiftyReg")

The mmand package for image processing may also be useful, and is used in some of the examples below.

Contents

Reading and writing images

RNiftyReg may be used to register and manipulate two and three dimensional images of any sort, although its origins are in medical imaging. Medical images in the standard NIfTI-1 format may be read into R using the readNifti function, which is based on the first-party RNifti package.

library(RNiftyReg)
image <- readNifti(system.file("extdata", "epi_t2.nii.gz", package="RNiftyReg"))

This image is an R array with some additional attributes containing information such as its dimensions and the size of its pixels (or voxels, in this case, since it is a 3D image).

As mentioned above, however, images do not have to be in NIfTI-1 format. Any numeric matrix or array can be used, and standard image formats such as JPEG and PNG can be read in using additional packages. For example,

## install.packages("jpeg")
library(jpeg)
image <- readJPEG(system.file("extdata", "house_colour_large.jpg", package="RNiftyReg"))

Complementary writeNifti and writeJPEG functions are provided by the RNifti and jpeg packages, respectively.

Image registration

Once two or more images have been read into R, they can be registered. Registration is the dual operation of searching a space of transformations for the best way to align two images, and then resampling one image onto the grid of the other. The images need not be the same size or in the same orientation.

There are two main classes of transformation available: linear and nonlinear. Linear, and specifically affine, transforms can represent translation, scaling and rotation in 2D or 3D space. They have up to 12 degrees of freedom and are appropriate to capture global shifts between images. Nonlinear transformations have many more degrees of freedom, and can capture localised distortions between images, but they are more time-consuming to estimate and more complex to work with.

Some sample 3D medical images are included with the package. We begin by registering two brain scan images, of the same person, with different contrasts. First we read them in, and then we pass them to the package's core registration function, niftyreg.

source <- readNifti(system.file("extdata", "epi_t2.nii.gz", package="RNiftyReg"))
target <- readNifti(system.file("extdata", "flash_t1.nii.gz", package="RNiftyReg"))
 
result <- niftyreg(source, target)

The last command will take a few seconds to complete. The result is a list with a number of components, the most important of which is the resampled source image in target space, result$image.

By default the transformation will be an affine matrix. If we want to allow for a nonlinear transformation, with scope for local deformations between one space and the other, we can perform an additional registration as follows.

result <- niftyreg(source, target, scope="nonlinear", init=forward(result))

Notice the scope argument, and also the fact that we use the result of the previous linear registration to initialise the nonlinear one. (The forward function extracts the forward transformation from the previous registration.) This should result in reduced convergence time, since the affine transformation provides a first approximation for the nonlinear registration algorithm. Nevertheless, this algorithm will generally take longer to complete.

Applying and manipulating transformations

Once a transformation between two images has been established through registration, it can be extracted, applied to another image or pixel coordinates, or manipulated. Transformations can also be read from or written to file, or created from scratch. Registration is by default symmetric, meaning that forward and reverse transformations are calculated simultaneously. These can be extracted using the forward and reverse functions.

Let's use a simple image by way of example. We will need the jpeg package to read it in, and the mmand package (version 1.2.0 or later) to visualise it.

## install.packages(c("jpeg","mmand"))
library(jpeg); library(mmand)
 
house <- readJPEG(system.file("extdata", "house_colour_large.jpg", package="RNiftyReg"))
display(house)

House photo

Clearly this is a colour image, with red, green and blue channels. RNiftyReg can work with it in this format, but internally the channels will be averaged before the registration starts. This step performs a colour-to-greyscale conversion equivalent to

house_bw <- apply(house, 1:2, mean)
display(house_bw)

Greyscale house photo

Now, instead of registering the image to another image, let's create a simple affine transformation that applies a skew to the image.

affine <- buildAffine(skews=0.1, source=house)
print(affine)
# NiftyReg affine matrix:
#  1.0  -0.1   0.0   0.0
#  0.0   1.0   0.0   0.0
#  0.0   0.0   1.0   0.0
#  0.0   0.0   0.0   1.0
# Source origin: (1, 1, 1)
# Target origin: (1, 1, 1)

So, this is a diagonal matrix with just a single off-diagonal element, which produces the skew effect. (The sign is negative because NiftyReg actually represents transforms from target to source space, not the more intuitive reverse.) Let's apply it to the image using the important applyTransform function, and see the effect.

house_skewed <- applyTransform(affine, house)
display(house_skewed)

Skewed house photo

Moreover, we can transform a pixel coordinate into the space of the skewed image:

applyTransform(affine, c(182,262,1))
# [1] 208.1 262.0   1.0

Notice that the skew changes the first coordinate (in the up-down direction), but not the second (in the left-right direction) or third (the colour channel number).

Finally, we can register the original image to the skewed one, to recover the transformation:

result <- niftyreg(house, house_skewed, scope="affine")
print(forward(result))
# NiftyReg affine matrix:
#  1.000845194  -0.099988878   0.000000000  -0.279240519
# -0.003345727   1.000585556   0.000000000   0.543151319
#  0.000000000   0.000000000   1.000000000   0.000000000
#  0.000000000   0.000000000   0.000000000   1.000000000
# Source origin: (1, 1, 1)
# Target origin: (1, 1, 1)

Notice that the estimated transformation closely approximates the generative one, with the element in row 1, column 2 being very close to -0.1. We can decompose this estimated transformation and recover the skew component:

decomposeAffine(forward(result))$skews
#        xy        xz        yz 
# 0.1032827 0.0000000 0.0000000 

Two other manipulations can be helpful to know about. The first is calculating a half-transform, which can be used to transform the image into a space halfway to the target. For example, using our registration result from above,

half_xfm <- halfTransform(forward(result))
display(applyTransform(half_xfm, house))

Half-skewed house photo

This results in half of the skew effect being applied. Finally, the composeTransforms function allows the effects of two transforms to be combined together. Combining a half-transform with itself will result in the original full transform.

all.equal(forward(result), composeTransforms(half_xfm,half_xfm), check.attributes=FALSE)
# TRUE

Convenience functions

The package provides a group of convenience functions—translate, rescale, skew and rotate—which can be used to quickly apply simple transformations to an image. For example, the skew operation applied above can be more compactly written as

house_skewed <- skew(house, 0.1)
display(house_skewed)

Skewed house photo

Since these take the image as their first argument, they are compatible with the chaining operator from the popular magrittr package. However, because such a chain applies multiple transformations to an image, there may be a loss of precision, or of data, compared to a single more complex operation. For example, while

library(magrittr)
house_transformed <- house %>% rotate(pi/4, anchor="centre") %>% translate(30)
display(house_transformed)

First transformed house

is much more readable than

xfm <- composeTransforms(buildAffine(angles=pi/4, anchor="centre", source=house), buildAffine(translation=30, source=house))
house_transformed <- applyTransform(xfm, house)
display(house_transformed)

Second transformed house

the latter avoids the creation of a black band across the top of the final image, since it has access to the full content of the original image, rather than just the truncated version produced by the rotation.

Upgrading to RNiftyReg 2.x

RNiftyReg 2.0.0 is a more-or-less complete rewrite of the package, with the goals of simplifying both the package's dependencies and its usage. The upstream NiftyReg code has also been updated. However, it should still be possible to read and use transformations created using RNiftyReg 1.x.

The core changes are

  • The oro.nifti package is no longer needed, nor used for reading and writing NIfTI files (RNiftyReg now offers readNifti and writeNifti, which are much faster). However, objects of S4 class nifti can still be used with the package if desired. Functions return either plain R arrays with attributes or bare-bones internalImage objects, which contain only some basic metadata and a pointer to a C-level data structure.
  • There are new functions for halving a transform (halfTransform), composing two transforms (composeTransforms), and building an affine transform from scratch (buildAffine).
  • Registration is now symmetric by default (for both linear and nonlinear), a newer symmetric nonlinear approach is now used, and default cost function weights have been tweaked. Therefore, the arguments to the core niftyreg function, and its linear and nonlinear special cases, have changed in both name and defaults. See ?niftyreg and related help pages for details.
  • It is no longer necessary to use functions specific to transform type to perform many operations. For example, the work of the old applyAffine, applyControlPoints, transformWithAffine and transformWithControlPoints functions is done by the flexible new applyTransform function. The forward and reverse transforms can always be obtained from a registration using the new forward and reverse functions, no matter what their type is. However, some affine-only functions, such as decomposeAffine, retain their names.
  • The affineType attribute has gone, and convertAffine is no longer a user-visible function. All affine matrices are stored using the NiftyReg convention. FSL-FLIRT affines can still be read in, but they are converted to NiftyReg convention immediately. In addition, source and target image information is attached to the transforms in attributes, and so does not need to be specified in most function calls.

News

Significant changes to the RNiftyReg package are laid out below for each release.

=================================================================================

VERSION 2.6.4

  • The algorithm for searching a deformation field is now more careful to avoid a potentially infinite loop.

=================================================================================

VERSION 2.6.3

  • The package now checks at load time whether it was built against a version of the RNifti package which is compatible with the currently installed one, and raises a warning prompting reinstallation if not.

=================================================================================

VERSION 2.6.2

  • Searching a deformation field is now done with gradient descent, rather than exhaustively. This makes applyTransform() substantially faster when applying a nonlinear transform to points.

=================================================================================

VERSION 2.6.1

  • The target image for a half transform now has its xform modified to better represent a half-way space between the original source and target.
  • The package now uses a configure script to more robustly detect OpenMP support.

=================================================================================

VERSION 2.6.0

  • New convenience functions translate(), rescale(), skew() and rotate() can be used to apply simple transformations to images easily.
  • The niftyreg() group of functions gain a "precision" argument, which controls the working precision of the optimisation.
  • A single angle supplied to buildAffine() with a 2D source file will now be interpreted in the way that would usually be expected (i.e., equivalent to a yaw angle).
  • Compositions of two affine transformations appear to have been happening in the wrong order. This has been corrected.
  • ANALYZE-format files, whose orientation is intrinsically ambiguous, are now assumed to be use LAS storage order (with a warning).
  • The registration functions will no longer modify the low-level C structures attributed to their arguments.
  • The "loder" package is now used in preference to "png".

=================================================================================

VERSION 2.5.0

  • RNiftyReg can now register 2D multichannel (RGB or RGBA) images. The red, green and blue channels are averaged internally before applying the algorithm. Previously such images would produce an error or segfault. (Reported by John Muschelli, issue #4.)
  • The decomposeAffine() function now uses a Cholesky decomposition internally. (Hat-tip to Tim Tierney.)
  • The "nBins" argument to niftyreg.nonlinear() was not previously passed to the NiftyReg library. This has been corrected.
  • Applying a scale transform to a new image with applyTransform() should no longer produce a warning. (Reported by Tim Tierney.)
  • The NiftyReg library has been updated.

=================================================================================

VERSION 2.4.1

  • The package has been updated for compatibility with "RNifti" version 0.3.0.

=================================================================================

VERSION 2.4.0

  • The package's functions for working with NIfTI-1 files and associated image data structures have been moved into the new "RNifti" package, since they are useful more broadly. RNiftyReg is therefore more focussed on registration.
  • The "internal" argument to niftyreg(), and its linear and nonlinear variants, can now be set to avoid creating any internal images in the result. This should only be necessary in exceptional circumstances, and the default behaviour remains as before.

=================================================================================

VERSION 2.3.0

  • The new similarity() function can be used to calculate a metric of similarity between two images (currently only normalised mutual information). This should not decrease after registration!
  • The composeTransforms() function now accepts more than two arguments, and handles NULL arguments (corresponding to no transform).
  • The package can now work with "MriImage" objects (from package "tractor.base").
  • It is now possible to interrupt the main registration functions.
  • The pixdim() and pixunits() functions are now (S3) generic.
  • The C++ code underlying the updateNifti() function could sometimes produce spurious errors about missing strings, when passed a list. This has been corrected. (Reported by Tim Tierney.)

=================================================================================

VERSION 2.2.0

  • The new asAffine() function can be used to coerce a numeric matrix to the "affine" class.
  • The updateNifti() function now accepts a named list of header elements, of the sort generated by dumpNifti(), as its second argument. This allows elements of the image's internal metadata to be updated manually, for advanced users.
  • It is now possible to disable the comment lines usually produced by writeAffine().
  • The niftyreg() function, and the two more specialised functions, gain an "internal" argument, which results in an internal image as output.
  • There is now a replacement method for pixunits().
  • The result of calling dumpNifti() now has a class and a print() method.
  • Elements of the "scales" argument to buildAffine() may now be negative.
  • Setting xforms would not work if the image did not already contain a valid internal pointer attribute. This has been corrected.
  • An incompatibility with C++11 has been fixed.

=================================================================================

VERSION 2.1.0

  • New functions for replacing an image's sform and qform matrices have been added. These are for advanced use only.
  • Some tests are skipped on Solaris, to avoid failures which have been hard to reproduce.

=================================================================================

VERSION 2.0.1

  • Fixes for Solaris compilation and UBSAN issues.

=================================================================================

VERSION 2.0.0

  • This is a more-or-less complete rewrite of the package, with the goals of simplifying both the package's dependencies and its usage. The upstream NiftyReg code has also been updated. However, it should still be possible to read and use transformations created using RNiftyReg 1.x.
  • The oro.nifti package is no longer needed, nor used for reading and writing NIfTI files (RNiftyReg now offers readNifti() and writeNifti(), which are much faster). However, objects of S4 class "nifti" can still be used with the package if desired. Functions return either plain R arrays with attributes, or bare-bones "internalImage" objects, which contain only some basic metadata and a pointer to a C-level data structure.
  • There are new functions for halving a transform (halfTransform), composing two transforms (composeTransforms), and building an affine transform from scratch (buildAffine).
  • Registration is now symmetric by default (for both linear and nonlinear), a newer symmetric nonlinear approach is now used, and default cost function weights have been tweaked. Therefore, the arguments to the core niftyreg() function, and its linear and nonlinear special cases, have changed in both name and defaults. See ?niftyreg and related help pages for details.
  • It is no longer necessary to use functions specific to transform type to perform many operations. For example, the work of the old applyAffine(), applyControlPoints(), transformWithAffine() and transformWithControlPoints() functions is done by the flexible new applyTransform() function. The forward and reverse transforms can always be obtained from a registration using the new forward() and reverse() functions, no matter what their type is. However, some affine-only functions, such as decomposeAffine(), retain their names.
  • The "affineType" attribute has gone, and convertAffine() is no longer a user-visible function. All affine matrices are stored using the NiftyReg convention. FSL-FLIRT affines can still be read in, but they are converted to NiftyReg convention immediately. In addition, source and target image information is attached to the transforms in attributes, and so does not need to be specified in most function calls.
  • Many smaller changes. The newly updated README file should be used as the main reference point for usage, supplemented by the R-level function documentation.

=================================================================================

VERSION 1.1.4

  • The package was not consistent about how to treat images without NIfTI orientation information, and this could lead to erroneous results when transforming points between spaces using an affine matrix. This has been corrected. (Reported by Jiří Borovec.)

=================================================================================

VERSION 1.1.3

  • The C++ code is now more careful to avoid buffer overruns when text fields in a "nifti" object are too long.

=================================================================================

VERSION 1.1.2

  • Some sample data has been added to the package, along with some code examples.
  • Code tweaks have been made for C++11 compatibility.

=================================================================================

VERSION 1.1.1

  • Under certain circumstances the NiftyReg C++ code could create a control point image corresponding to an invalid "nifti" object, resulting in a spurious error when running niftyreg.nonlinear(). This has been corrected. (Reported by Andreas Heindl.)

=================================================================================

VERSION 1.1.0

  • A set of functions has been added for transforming points between spaces using an affine or nonlinear transformation.
  • All registration functions gain an "estimateOnly" option, which can be used to estimate a transformation without actually resampling the image.
  • The new getDeformationField() function can be used to obtain the deformation field and/or Jacobian determinant map corresponding to a transformation.
  • The invertAffine() function has been added to invert an affine matrix while preserving its "affineType" attribute.
  • RNiftyReg will now use OpenMP, if supported by the current build of R.
  • Image arguments are now coerced automatically to the "nifti" class if they are not already of this type.
  • Image data is no longer copied between data structure types when not needed, saving a bit of time and memory.
  • It is no longer an error in readAffine() if the affine type is not specified or stored.

=================================================================================

VERSION 1.0.2

  • A further compiler-sensitive problem, which was causing the package to fail to install using Solaris Studio, has been addressed.

=================================================================================

VERSION 1.0.1

  • A Windows-specific compile-time problem has been addressed.

=================================================================================

VERSION 1.0.0

  • RNiftyReg now provides an interface to the symmetric version of the NiftyReg nonlinear registration algorithm. This allows consistent transformations to be obtained in both directions, in one operation.
  • The upstream code has been updated to match NiftyReg v1.3.9. Default cost function weightings have also been tweaked to match those used in that release.
  • The type of an affine matrix is now written to file with it. The type therefore does not necessarily need to be specified when reading the matrix back in.
  • The "interpolationPrecision" option has been removed where it appeared. The code now simply makes an appropriate choice in each case.
  • Internal calculation of the final control point spacing from an initial control point image was sometimes wrong. This has been corrected.

=================================================================================

VERSION 0.6.2

  • Nonlinear registration now works properly for 2D target images.
  • Control point spacing is now set from the initial control point image (as it should be), if one is supplied.
  • Two small bugs in the upstream NiftyReg code have been fixed.
  • Attempting to perform 4D-to-2D registration now produces a more specific error message.

=================================================================================

VERSION 0.6.1

  • The package now redefines standard C print functions to avoid problems with output not being captured by R.

=================================================================================

VERSION 0.6.0

  • RNiftyReg can now perform nonlinear registration via the niftyreg.nonlinear() function. This method optimises over a set of moveable control points arranged in a grid over the image. The applyControlPoints() function has also been added as the nonlinear equivalent of applyAffine(). See ?niftyreg.nonlinear for details on these new functions. The niftyreg() function now serves as a common interface to both the linear and nonlinear optimisations.
  • The "initAffine" argument to niftyreg() now accepts a list of matrices as well as a single affine matrix.
  • The upstream NiftyReg code has been updated, thereby taking in a number of tweaks added since the last release version.
  • A bug in the quaternion-to-matrix conversion in xformToAffine() has been fixed.

=================================================================================

VERSION 0.5.0

  • The niftyreg() and applyAffine() functions gain an "interpolationPrecision" argument, which affects the data type of the final image. See ?niftyreg for details.
  • The actual number of algorithm iterations run at each level is now stored in the output produced by niftyreg().
  • A bug has been addressed whereby niftyreg() returned the affine matrix used for initialisation, if specified, rather than the final one. (Reported by Brandon Whitcher.)
  • The first level of the algorithm is supposed to include both a rigid-body and affine optimisation, but the latter was previously omitted due to a bug. This has been corrected. In addition, the first level is allowed twice as many iterations as subsequent levels, as in the original NiftyReg implementation.
  • The "data_type" slot of the "nifti" class is now properly set in the image returned by niftyreg().
  • The underlying NIfTI reference library files have been updated.
  • Spurious parameters have been removed to avoid warnings at compile time.
  • Calls to problematic C/C++ functions such as exit() have been avoided (although NIfTI and NiftyReg library functions still use printf() and the like).

=================================================================================

VERSION 0.4.1

  • The reportr package is now properly imported in the namespace file.

=================================================================================

VERSION 0.4.0

  • The applyAffine() function has been added as a convenience wrapper to niftyreg(...,nLevels=0), for applying a precomputed affine transformation to a new image without further optimisation.
  • RNiftyReg now uses the reportr package for message and error reporting.
  • The package now has a formal namespace.

=================================================================================

VERSION 0.3.1

  • A bug in previous releases affecting the downsampling of images during early "levels" of alignment has been fixed. Coregistration of larger images, in particular, should be improved as a result.

=================================================================================

VERSION 0.3.0

  • 2D-to-2D, 3D-to-2D and 4D-to-3D registration are now supported by niftyreg(). As a result of this added flexibility, a list of affine matrices, rather than a single matrix, is now returned by this function.
  • An image of fewer than 4 pixels/voxels in any dimension cannot be registered since it cannot accommodate a single matching block. Providing such an image as the source or target now produces an error rather than a crash. (Thanks to Takeo Katsuki for pointing this out.)
  • Further documentation improvements.

=================================================================================

VERSION 0.2.0

  • Added an option to niftyreg() to control the type of interpolation that is applied to the final coregistered image.
  • Performing better type checking of image data before passing it to the C++ code, to avoid intermittent errors from R.
  • Documentation improvements.

=================================================================================

VERSION 0.1.0

  • First public release. 3D affine and rigid-body registration are available.

=================================================================================

Reference manual

It appears you don't have a PDF plugin for this browser. You can click here to download the reference manual.

install.packages("RNiftyReg")

2.6.5 by Jon Clayden, 16 days ago


https://github.com/jonclayden/RNiftyReg


Report a bug at https://github.com/jonclayden/RNiftyReg/issues


Browse source code at https://github.com/cran/RNiftyReg


Authors: Jon Clayden [cre, aut] , Marc Modat [aut] , Benoit Presles [aut] , Thanasis Anthopoulos [aut] , Pankaj Daga [aut]


Documentation:   PDF Manual  


Task views: Medical Image Analysis


GPL-2 license


Imports Rcpp, RNifti, ore

Suggests jpeg, loder, mmand, testthat

Linking to Rcpp, RcppEigen, RNifti


Imported by FIACH, patternize.

Suggested by spant.


See at CRAN