A 3D geometric transformation for a nonrigid image registration method
Norma Pilar Castellanos-Abrego*
* Universidad Autónoma Metropolitana -Iztapalapa, Electrical Engineering Department, D.F., 09340, México.
Correspondence: ]]>
Norma Pilar Castellanos Abrego.
xanum.uam.mx
Received article: 25/may/2010.
Accepted article: 30/november/2010.
ABSTRACT
A 3D geometric transformation is introduced for the nonrigid registration of medical images as an extension of a previous work carried out for two dimensions. A 3D spatial transformation is analyzed in order to guaranty the continuity, the differentiability and the one-to-one transformation by imposing constraints to the transformation parameters. It is also shown and analyzed the results when the fully automatic nonrigid registration method is applied to a CT-PET stack of the thorax with a spatial resolution of 80 x 80 x 80 and to a RM head stack with a spatial resolution of 128 x 128 x 128 pixels. The 3D geometric transformation has a spherical domain and it allows the continuity of the transformation in its boundary. This geometrical transformation can be applied to global or local ROIs (region of interest) up to a minimum diameter of three pixels. The nonrigid image registration method employs an evolutionary algorithm to obtain satisfactory global solutions while it maximizes the normalized mutual information (NMI). This approach has the disadvantage that the speed of convergence and the accuracy of the method depend on the population size of the evolutionary algorithm. Results show an improvement in the global similarity function between the target and source volumes throughout 73 transformations, from coarse to fine (3 levels of resolution), from 0.501 7 to 0.5033, using a population size of 10 individuals. 3D surface reconstructions of the thorax are also shown before and after the nonrigid registration. In addition, a simulated experiment is carried out with a RM head stack, where a unique transformation was applied. Here, it was got an improvement in the similarity criterion from 0.5046 to 0.5218.
Key words: Nonrigid image registration, nonlinear geometrical transformation.
]]> RESUMEN
Una tranformación geométrica tridimensional (3D) se presenta para el registro de imágenes médicas como la extensión de un trabajo previo desarrollado para dos dimensiones. La tranformación espacial 3D es analizada con la intención de garantizar la continuidad, la diferenciabilidad y la relación uniquívoca mediante la imposición de restricciones a los parámetros de la tranformación. Se demuestran y analizan los resultados cuando el registro no rígido completamente automático es aplicado a un conjunto de imágenes de tórax con una resolución espacial de 128x128x 128 pixeles adquiridas por el método de la tomografía computarizada asociada a la tomografía de emisión de positrones (PET-CT). La tranformación geométrica 3D tiene un dominio esférico que permite la continuidad de la transformación en su frontera. Esta transformación geométrica puede ser aplicada a regiones de interés locales o globales hasta un diámetro mínimo de 3 pixeles. El método de registro no rígido utiliza un algoritmo evolutivo para obtener resultados globales satisfactorios mientras que maximiza la información mutua normalizada. Esta propuesta tiene la desventaja de que la velocidad de convergencia y la precisión del método depende del tamaño de la población en el algoritmo evolutivo. Los resultados demuestran una mejora en la función de similitud global entre los volúmenes fuente y destino a través de 73 transformaciones, desde burdos hasta finos (en 3 niveles de resolución), de 0.5017 a 0.5033 usando un tamaño de población de 10 individuos. Las reconstrucciones 3D del tórax también se muestran antes y después de la transformación no rígida. Además, experimentos de simulación fueron desarrollados con imágenes tomadas de un estudio del cráneo por resonancia magnética utilizando sólo una transformación. Aquí, se obtuvo una mejora en el criterio de similitud desde 0.5046 a 0.5218.
Palabras clave: Registro de imagen no rígido, transformación geométrica no lineal.
INTRODUCTION
Mathematical models for the nonrigid image registration are divided into two main categories: the physically based models (e.g., linear elasticity and fluid flow); and the basis function expansions (e.g., radial basis functions, B-splines and wavelets). If information about geometric differences between images is available, that information should be used to select the transformation. If no such information is available, a transformation function that can adapt to local geometric difference between the images should be chosen. A vast review of geometric transformations is carried out by L. Zagorchev in 20061 and for M. Holden in 20081 where they report the advantages and disadvantages of the different approaches of nonrigid medical image registration procedures. They compare the performance of the nonrigid image registration methods, finding differences based on the size of the set of control points and the length of spacing between them.
Recently there have been proposed other methods for the nonrigid image registration. One of them demonstrates advantages in reducing the degrees of freedom of the transformation without losing accuracy, nevertheless, it strongly depends on the application3. Thus, research goes on working on validation methods for particular applications: with brain EPI-MRI, authors compare a number of similarity functions as well as the statistical results obtained from brain anatomy expert's evaluation4, the semiautomatic registration of pre-and postbrain tumor resection5 and the intensity standardization of MRI6.
Other works are related to accuracy and automatism6, the speed of convergence6 or robustness7.
In this work a 3D geometrical transformation is introduced for the nonrigid registration of medical images. As an extension of the previous work reported in4 in order to guaranty the continuity, the differentiability and the one-to-one transformation by imposing new constraints to the transformation parameters. Besides, we applied the transformation to a volume of medical images in order to assess the performance of the transformation as part of a nonrigid registration algorithm.
This approach is motivated for the simplicity of the transformation in the sense that it employs a reduced searching space in comparison with other methods where its dimensionality is increased as locality of the transformation is required (reduced space between control points). In our case, a spherical domain or region of interest (ROI) is placed where the local or global transformation will be applied; it has the advantage over other methods that the optimization algorithm only searches for 10 parameters. It is useful for correcting small deformations as for example in the image fusion of SPECT-CT where nonrigid registration is needed to alleviate the movements during the image acquisition or due the breathing, our next goal in this research. For now, we are introducing the 3D transformation, the analyses carried out and the optimization method employed. Besides, the transformation allows the application of local and global transformations by the composition of functions. On the other hand, the performance of the similarity function is tested, the normalized mutual information (NMI) that has been extensively used, demonstrated its accuracy and robustness in rigid image registration, its extension to nonrigid image registration is not trivial and it has been reported as an active field of research8.
]]>METHODOLOGY
Let Ω and Ω, be the source and target volumes domains, respectively, where Xs ∈ Ωs and Xt∈ Ωt.
The proposed 3D geometrical transformation g: ℜ3 → ℜ3 (eq. 1) carries out a smooth map from the center of a sphere to its boundary, by applying an affine transformation A , and the identity I , at the center and at the boundary, respectively. In this way the continuity of the transformation is guarantied at the boundary. α∈ℜ is called the smoothing parameter because it establishes the differentiability of eq. 1, the 9 matrix coefficients in A are real numbers (eq. 2).
Function λ (eq. 3), is a norm vector ||•|| in ℜ3 which defines the spherical domain, where Xes and Xmax, are the central position and the boundary of the sphere, respectively. The spherical domain is normalized such that. ||Xes -Xs|| ≤ 0.5 From eq. 3, it can be analyzed that there must be at least one pixel between the center and the boundary, resulting in a transformation domain with a minimum diameter of 3 pixels. Then, eq. 3 can be written as λ(Xs) = 2((0.5-Xs)2 + (0.5 - ys)2 + (0.5 - zs)2)1/2- In the other hand, the smoothness parameter must be constrained to be greater than 2 to avoid singularities in the center of the sphere.4 Because of in the boundary the identity transformation is applied, the continuity here is of order zero (C0).
The behavior of the 3D geometrical transformation can be observed by means of the intersection of the contour surfaces of eq. (1) into the spherical domain. Figure 1 shows the transformation of the center of the normalized sphere (xs = (0.5, 0.5, 0.5)). In this case, Xs is lightly moved to a new position that is defined by the vector of parameters (p). This vector of parameters establishes a one-to-one correspondence from the positions in the source volume to the transformed ones. It means that after the application of the geometrical transformation, a source position is moved to a new and unique target position. This property is not always achieved by the geometrical transformation as it is shown below.
In figures 2, 3 and 4 are shown violations to the one-toone correspondence described by the vector of parameters p. It can be observed in these figures how the contour surfaces intersect each other more than once or an intersection is absent. For example, in figures 2 and 4 the position in the source volume Xs = (0.95, 0.5, 0.8), is transformed to several new positions. This is due to the surface contour curvature always defined by the vector of parameters. In figure 4 is shown an example of a lack of a junction point among the three contour surfaces in the spherical domain.
]]>
One-to-one transformation. As was mentioned and shown before, some constraints must be imposed in order to guaranty a one-to-one correspondence and therefore a valid geometrical transformation. It is carried out by analyzing the the maximum and minimum slopes by planes of the countour surfaces in such a way that each surface curvature will be limited (eqs. (4) to (6)).
where
The non-rigid image registration algorithm. The 10 parameters of the proposed 3D geometrical transformation are optimized by an evolutionary algorithm that searches for a global solution. It maximizes the normalized mutual information (NMI) and incorporates the constraints described in eqs. (4) to (6).
The region of interest is selected automatically, the first transformation is applied to a global ROI and then the volume is divided each i time by two generating 8 new local volumes. This process of localizing local ROIs in each i level goes on until a wished diameter is reached. The composition of functions is applied throughout each transformation and a trilinear interpolation is used.
Optimization algorithm. A hybrid genetic algorithm (HGA) is used to search for the optimal global parameters, p. The hybrid algorithm reported here is an adaptation of that reported by Pham et al9, which has a variable mutation rate.
Let fs and ft, be the intensity distribution of the source and target images, respectively. The optimization algorithm is formulated in the following way: Maximize the NMI,
]]>
subject to a search space S∈ℜ10 and a feasible region F ⊆ S, defined by the 3 constraints described in eqs. (4) to (6) and α > 2.
Let Q (r) be a population of N individuals uniformly distributed and conformed by the parameters of the transformation, 
in the r-th generation, which are initialized to be in a feasible region (F).
1. The /-th value of the objective function 
is assigned, where the best vector of parameters pb (r) with the best cost function Ib (r) is selected.
2. The reproduction of individuals is done by assigning p/ to the next generation as described as follows,
where the relation 
denotes a positive weight to improve the performance of the HGA. It depends on the average value of p. The reproduction is not carried out if the new generation lies outside the feasible region.
3. Mutation is regulated by the mutation rate operator for the r generation, choosing randomly one of the parameters, with normal probability function. The mutation is carried out with the expression
where var is a random variable and u the parameter to be mutated.
]]> 4. The HGA ends when a number of iterations with the same solution is reached (rep_max). Nevertheless, a relationship is established between each c iteration, the initial and maximum mutation rate: mutrate_min and mutrate_max, respectively, in order to determine the current mutation rate as follows:mutrate = (c-1)* mutrate _max/rep _max + mutrate_min.
The constant values are found empirically: mutrate_min = 0.01, mutrate_max = 0.5 , N = 10 and rep_max = 10.
RESULTS AND DISCUSSION
The whole nonrigid image registration algorithm was implemented in MATLAB as well as the 3D renderings. This algorithm has been applied to a CT and PET stacks of axial images of the thorax with a resolution of 80 x 80 x 80 pixels. These stacks were put in the center of the spherical domain generating a manual rigid registration. The CT stack was the reference or target volume (Figura 5), and the PET transmission volume is the source volume (Figura 6 (a)). After the parameters are optimized, they are applied to the PET emission volume.
The algorithm was able to improve the global similarity (NMI) between the target and source volumes from 0.5017 to 0.5028 for the PET transmission volume, in the first level (Figura 7) and to 0.5033 throughout 3 levels or 73 transformations. It was found that a good rigid registration is essential in order to obtain acceptable results. Also, for a global transformation, it must be guaranteed that all the volume is inside the spherical domain. In order to show the validity of the geometrical transformation of the above example, the optimized global vector of parameters is represented in the countour surfaces shown in Figure 8. It can be observed from the curvature of the three contour surfaces and their intersection, that there is a unique transformed position for the source point, indicating a valid geometrical transformation.
]]>
Our last experiment is for a head of MR sagittal images (Figura 9 (a)) with a simulated deformation (Figura 9 (b)) with the vector of parameter p=[-0.34, 0.1, 0.2, 0.2, -0.2, -0.03, 0.16, 0.02, -0.46, 2.15]. Here, the method was able to improve the NMI from 0.5046 to 0.5218 choosing the volume with simulated deformation as the target. In Figura 10 is shown the resulting 3D rendering which is very similar (qualitatively) with the expected volume (Figura 9 (b)). In fact, the resulting vector of parameters was p=[-0.19, 0.06, 0.09, -0.04, -0.3, 0.19, 0.07, 0.05, -0.42, 2.08], demonstrating the difficulty of finding the maximum of the NMI in a very close region into the searching space. Nevertheless, the evolutionary algorithm was able to find an acceptable solution.
CONCLUSIONS
A 3D geometrical transformation is introduced in this work as part of a nonrigid image registration algorithm. This algorithm was able to improve the global similarity criterion by applying the composition of functions and optimizing the parameters of each transformation. As an advantage of the proposed geometric transformation, it generates a low dimension searching space diminishing the complexity of the numerical solution. Nevertheless, the transformation has the disadvantage that is not invertible.
There has been analyzed the constraints of the transformation in order to guaranty its applicability as a nonrigid registration of medical volumes. We have demonstrated its adequate performance in two volumes of different image modality. More tests must be done with other image modalities to prove its robustness.
For accuracy comparison with other methods, there must be done more testes but focusing in a particular application, with standard data and similarity criterions, and with the assessment of an expert physician, which is not the objective at the present. Nevertheless, we demonstrate in this work, that the method for 3D nonrigid image registration improves the local and global similarity function throughout different levels of locality.
Although, it was able to improve the NMI, future work will be related with the incorporation of other kind of similarity functions as the conditional mutual information reported in.8
]]>REFERENCES
1. Zagorchev L, Goshtasby A. A Comparative Study for Nonrigid Image Registration", IEEE Trans. Image Processing 2006; 15(3): 529-538. [ Links ]
2. Holden M. A Review of Geometric Transformations for Nonrigid Body Registration. IEEE Trans Medical Imaging 2008; 27(1): 111-128. [ Links ]
3. Nanayakkara ND, Chiu B, Samani A, Spence JD, Samara-bandu AF. A twisting and bending Model-Based Nonrigid Image registration Technique for 3-D Ultrasound carotid Images. IEEE Trans Medical Imaging 2008; 27(10): 13781388. [ Links ]
4. Gholipour AN, Kehtarnavaz RW, Briggs KS, Gopinath W, Ringe A, Whittemore S, Cheshkov K. Bakhadirov. Trans Biom Eng 2008; 55(2): 563-571. [ Links ]
]]>5. Ding S, Miga MI, Noble JH, Cao A, Dumpuri P, Thompson RC, Dawant BM. Semiautomatic Registration of Pre-and Post-rain Tumor Resection Laser Range Data: Method and Validation. Trans Biom Eng 2009; 56(3): 770-780. [ Links ]
6. Jäger F, Hornegger J. Nonrigid Registration of Joint Histograms for Intensity Standartization in Magnetic Resonance Imaging. IEEE Trans Medical Imaging 2009; 28(1): 137-150. [ Links ]
7. Sdika M. A Fast Nonrigid Image Registration with Constraints on the Jacobian using Large Scale Constrained Optimization. IEEE Trans Medical Imaging 2008; 27(2): 271-281. [ Links ]
8. Loeckx D, Slagmolen P, Vandermeulen D, Suetens P. Nonrigid Image Registration Using Conditional Mutual Information. Trans Medical Imaging 2010; 29(1): 19-29. [ Links ]
9. Pham DT, Karaboga D. Intelligent Optimization Techniques Springer-Verlag, Londres 2000: 51-61. [ Links ]
]]>Nota
Este artículo también puede ser consultado en versión completa en: http://www.medigraphic.com/ingenieriabiomedica/
]]>
