About this Journal Submit a Manuscript Table of Contents
Computational Intelligence and Neuroscience
Volume 2012 (2012), Article ID 561406, 7 pages
http://dx.doi.org/10.1155/2012/561406
Research Article

Hybrid Particle Swarm Optimization and Its Application to Multimodal 3D Medical Image Registration

College of Information and Science, Ritsumeikan University, Shiga 525-8577, Japan

Received 27 April 2012; Accepted 3 June 2012

Academic Editor: Huiyan Jiang

Copyright © 2012 Chen-Lun Lin et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

In the area of medical image analysis, 3D multimodality image registration is an important issue. In the processing of registration, an optimization approach has been applied to estimate the transformation of the reference image and target image. Some local optimization techniques are frequently used, such as the gradient descent method. However, these methods need a good initial value in order to avoid the local resolution. In this paper, we present a new improved global optimization approach named hybrid particle swarm optimization (HPSO) for medical image registration, which includes two concepts of genetic algorithms—subpopulation and crossover.

1. Introduction

In the area of medical image analysis, multimodality 3D image registration is an important issue [1]. The purpose of image registration is to register a target image (moving image) to a reference image (fixed image) so that we can combine the information of two images to obtain more detailed information or some specific features. For example, the PET image usually shows metabolic activity of organs and abnormal tissues clearly but lacks the texture of organ tissues. On the other hand, the MR image is described by much complex intensity to represent the texture of organ tissues well. If we implement the MR-PET image registration to combine the information of two images which are different modality, then we can get the accurate shape, volume, and location of abnormal tissues from the registered image. The registration is a very important and helpful preprocessing technique for medical diagnosis or surgical operations.

The processing of registration can be seen as an iterated optimization framework, and it can be divided into 3 parts: transformation, cost function, and optimization. In each iteration, the target image is firstly transformed by transformation according to the parameter of the current time. Then, the reference image and the transformed target image are used to calculate the cost function which can evaluate whether the two images are registered or not under the current parameter of transformation. If the images are not registered, an optimization method will be used to adjust the parameters, and a new iteration will start.

The application of registration can be classified with dimensionality of image, modality of image, and model of transformation. There are 2D to 2D, 3D to 2D, and 3D to 3D image registration for many different application. The 3-D to 3-D image registration usually needs to estimate more parameters than the 2-D to 2-D image registration, so it does require a more advanced optimization method. Then, if two images which have different scope of intensity have to be registered, such as CT-MR registration, it is called multimodal registration. On the other hand, if register two images which have same modality, it is called monomodal registration, such as CT-CT. Depending on the modal of registration, different cost functions have been used, such as the sum of squared intensity difference (SSD) for mono-modal registration or mutual information (MI) for multimodal registration. Moreover, the type of transformation model determines whether a registration belongs to rigid or nonrigid one. If target objects which we want to register are different in shape or deformable such as liver, the nonrigid transformation is used it is called nonrigid registration [2]; otherwise it is called rigid registration. This paper is focused on rigid multimodal 3D medical image registration.

As our previous explanation, we estimate the parameter of transformation by optimizing a cost function (similarity metric) in the processing of registrations. Some local optimization techniques, such as the gradient descent method, are frequently used for medical image registrations [3, 4]. However, since the transformation parameters are generally nonconvex and irregular, these kind of methods require very good initial values in order to avoid the local minimum.

To overcome the local resolution problem, the genetic algorithm (GA), which is one of the global optimization techniques, has been proposed for medical image registrations [5]. Although GA is an advanced method for global optimization, it requires huge computation time and lacks the fine tuning capabilities. We need more powerful approach. Particle swarm optimization (PSO) is a new global optimization technique. This method is a stochastic, population-based evolutionary computer algorithm [6, 7]. PSO is an extremely simple algorithm, and it seems to be more effective for optimizing a wide range of functions, and has been shown very effective for 2D rigid image registration [8].

On the other hand, more transform parameters must be estimated in 3-D image registrations. Our experiments show that conventional GA and PSO cannot find the global optimum well. Thus, we propose a new approach named hybrid particle swarm optimization (HPSO) for PSO. In our proposed method, two concepts of genetic algorithms—subpopulation and crossover—are incorporated into the conventional PSO method to improve the accuracy of that conventional GA and PSO, because these conventional methods can not find the global optimum resolution when we need estimated a huge number of parameters. Experiments are done with both mathematical test functions and medical volume data, and it is demonstrated that the proposed HPSO performs much better results than conventional gradient decent, GA, and PSO methods.

The paper is organized as follows: the image registration technique is summarized in Section 2, the particle swarm optimization (PSO) and our proposed hybrid particle swarm optimization (HPSO) are presented in Section 3, the experimental results with both mathematical test functions and medical volume data are presented in Section 4, and finally the conclusion and future works are given in Section 5.

2. 3D Image Registration

Medical Image registration is one of the fundamental tasks within medical image processing. The framework of medical image registration is shown in Figure 1. Two volumes of medical data that have to be registered are given as the fixed image and the moving image. The fixed image is denoted by 𝑓1(𝐱), where 𝐱 is a set of coordinates. The moving image is similarly denoted by 𝑓2(𝐱). Given 𝐓 is a transformation from the coordinate frame of the fixed image to the moving image, 𝑓2(𝐓(𝐱)) is the moving image associated with fixed image 𝑓1(𝐱). In order to simplify some of the subsequent equations, we will use 𝐓 to denote both the transformation and its parameterization. Here, we used an estimation of the transformation that registers the fixed image and moving image by maximizing their cost function (similarity metric) as shown in (1): 𝐓=argmax𝐓𝑓Mtric1(𝐱),𝑓2,(𝐓(𝐱))(1) where 𝐱 is the coordinate of a 3D or 2D point.

561406.fig.001
Figure 1: Framework of medical image registration.
2.1. Transformation

Registration can be seen as a process of finding the spatial transformation that maps points from one image to the corresponding points in another image. Here, we focus on rigid 3-D global transformation. The rigid transformation deals with 6 degrees of freedom for 3-D object translation and rotation. The rigid transform for 3-D object can be expressed by (2):𝐓Global𝑥𝑦𝑧+𝑡(𝐱)=𝐑𝐱+𝐭=cos𝛽cos𝛾cos𝛼sin𝛾+sin𝛼sin𝛽cos𝛾sin𝛼sin𝛾cos𝛼sin𝛽cos𝛾cos𝛽sin𝛾cos𝛼cos𝛾sin𝛼sin𝛽sin𝛾sin𝛼cos𝛾cos𝛼sin𝛽sin𝛾sin𝛽sin𝛼cos𝛽cos𝛼cos𝛽𝑥𝑡𝑦𝑡𝑧,(2)

where α, β, γ are rotation angles around each axis and 𝑡𝑥, 𝑡𝑦, 𝑡𝑧 are translations around each axis, respectively. It should be noted that there are only two parameters to be estimated for 2D rigid transform.

2.2. Metric Function (Mutual Information)

For multimodality medical image registration, mutual information (MI) is widely used as a similarity metric [9]. MI is an intensity-based similarity measure and is closely related to joint entropy. Given an image 𝐴 and an image 𝐵, the joint entropy of two images can be calculated by 𝐻(𝐴,𝐵)=𝑎,𝑏𝑝𝐴,𝐵(𝑎,𝑏)log𝑝𝐴,𝐵(𝑎,𝑏),(3) where 𝑝𝐴,𝐵 is the joint probability distribution function of pixels associated with images 𝐴 and 𝐵. The joint entropy is minimized when there is a one-to-one mapping between the pixels in 𝐴 and their counterparts in 𝐵. It increases while the statistical relationship between 𝐴 and 𝐵 weakens. Mutual information can be defined in terms of entropy as follows: 𝑀𝐼(𝐴,𝐵)=𝐻(𝐴)+𝐻(𝐵)𝐻(𝐴,𝐵),(4) where 𝐻(𝐴) and 𝐻(𝐵) are the individual entropies, which can also be represented by the probability distribution function (PDF) as 𝐻(𝑋)=𝑥𝑝𝑋(𝑥)log𝑝𝑋(𝑥).(5) Usually, a discrete joint histogram is adopted to estimate the joint PDF for the calculation of MI. Attempting to find the most complex overlapping regions, we should maximize the individual entropies and minimize the joint entropy which could explain each other well.

3. Hybrid Particle Swarm Optimization

In this paper, we propose a new approach named hybrid particle swarm optimization (HPSO) for 3-D rigid medical volume registration. For describing more details, we show the traditional PSO first.

3.1. Particle Swarm Optimization (PSO)

Assume an extent diffuse population existing which is called a swarm and an individual member of the swarm is termed as a particle. A particle can be imaged as a point in search space. A group of particles tend to cluster at a position where optimized results are identified. Therefore, to achieve particle swarm optimization, each particle adjusts itself by comparing previous experience and its neighbors to obtain the best result. The formula of particle swarm optimization can be represented as below: 𝐯𝑖𝑡+1=𝑤𝑡𝐯𝑡𝑖+𝑐1rand𝐩_𝐛𝐞𝐬𝐭𝑖𝐱𝑡𝑖+𝑐2𝐠rand𝐛𝐞𝐬𝐭𝐱𝑡𝑖,𝑤𝑡+1=𝑤𝑡𝑤+𝑑𝑤;𝑑𝑤=min𝑤max𝑇,𝐱𝑖𝑡+1=𝐱𝑡𝑖+𝐯𝑖𝑡+1.(6)

At iteration 𝑡, 𝐱𝑖 is the 𝑖th particle that moves with a velocity vector 𝐯𝑖, 𝐩_𝐛𝐞𝐬𝐭𝑖 is the personal best of 𝐱𝑖, and 𝐠_𝐛𝐞𝐬𝐭 is the global best among all particles. 𝑤 means the weight. The initial weight is 0.4. The minimum weight is 0.4 and the maximum is 0.98. 𝑐1 and 𝑐2 represent acceleration constants (𝑐1=𝑐2=2.0). rand is a uniformly distributed random number among 0 to 1. The movement of each particle is shown in Figure 2.

561406.fig.002
Figure 2: The movement of each particle.

The flowchart of PSO for image registration is shown in Figure 3.

561406.fig.003
Figure 3: Flowchart of PSO.

In the PSO-based registration approach, the particle 𝐱 is the transform parameter that needs to be estimated, the p_best is the maximum MI (cost function) of each particle, and g_best is determined by the cluster MI. MI in Figure 2 represents the similarity metric function.

3.2. Hybrid Particle Swarm Optimization

In this section, we propose a new approach named hybrid particle swarm optimization (HPSO) for 3-D rigid medical volume registration. Our method incorporates two concepts of genetic algorithms, which are subpopulation and crossover, into the traditional PSO. We expect that our proposed method will improve the accuracy of registration by taking advantage of subpopulation and crossover. The flowchart of HPSO is shown in Figure 4.

561406.fig.004
Figure 4: A total of 40 runs of 100 generations each for each function and each method were performed. The population size is 56 for GA, PSO, and HPSO, respectively, and the number of subpopulation is 4 for HPSO. The relative mean square error between the estimated solutions and real solutions is calculated for each trial. Table 1 shows the averaged relative mean square error over 40 trials for each method. Flowchart of HPSO.
3.2.1. Subpopulation

The particles are divided into a number of subpopulations. Each subpopulation has its own best optimum, 𝑔sub-best𝑘. The process of PSO is done for each subpopulation group. If the 𝑔sub-best𝑘 is better than 𝑔_best, then 𝑔_best, is replaced by the 𝑔sub-best𝑘, where 𝑘 is the subpopulation number.

3.2.2. Crossover

The 𝑔sub-best𝑖 are sorted in order with large mutual information. The top two 𝑔sub-best are selected as parents (𝐱𝑖 and 𝐱𝑗) for crossover, where 𝑖 and 𝑗 are their subpopulation number. The offspring are generated for each by arithmetic crossover, which are shown as 𝐱𝑖=rand𝐱𝑖+(1rand)𝐱𝑗,𝐱𝑗=rand𝐱𝑗+(1rand)𝐱𝑖,(7) and the velocities are given by 𝐯𝑖=𝐯𝑖𝐯𝑉,𝑗=𝐯𝑗𝐯𝑉,𝑉=𝑖+𝐯𝑗𝐯𝑖+𝐯𝑗,(8) where rand is a uniformly distributed random number among 0 to 1. The worst particle in the same subpopulation is replaced by the offspring.

4. Experimental Results

In this section, we perform several registration experiments with both test functions and medical volume data to evaluate the performance of the proposed HPSO technique. In the meantime, we also perform conventional GA and PSO for comparison. Finally, we implement the parallel technique to reduce computation cost and show the efficiency of our implementation.

4.1. Test Function Evaluation

In the first part of Section 4, we apply 4 mathematical functions which have different numbers of parameters to estimate the accuracy of our proposed method. The functions which were used are shown as below:𝑥𝐹1𝑓𝑖𝑖=1,,3=3𝑖=1𝑥2𝑖,𝑥𝑖[],𝑥5.12,5.11(9)𝐹2𝑓𝑖𝑖=1,2𝑥=10021𝑥2+1𝑥212,𝑥𝑖[],𝑥2.048,2.047(10)𝐹3𝑓𝑖𝑖=1,,20=(20×10)+20𝑖=1𝑥2𝑖10cos2𝜋𝑥𝑖,𝑥𝑖[],𝑥5.12,5.11(11)𝐹4𝑓𝑖𝑖=1,,10=1+10𝑖=1𝑥𝑖2400010𝑖=1𝑥cos𝑖𝑖,𝑥𝑖[].5.12,5.12(12)

The functions are commonly used for evaluation experiments and cover a variety of characteristics that affect algorithmic performance [10]. 𝐹1 is unimodal with global minimum at center and has 3 parameters to be estimated. 𝐹2 is strong epitasis with global minimum at center and has 2 parameters to be estimated. 𝐹3 is highly multimodal with global minimum at center and has 20 parameters to be estimated. 𝐹4 is multi-modal with global minimum at corner and has 10 parameters to be estimated. 𝐹3 and 𝐹4 are more difficult than 𝐹1 and 𝐹2, and 𝐹3 is the most difficult problem.

As shown in Table 1, both PSO and HPSO can get perfect solutions, and GA can also get a reasonable result for 𝐹1 and 𝐹2 in which the number of parameters is only 3 and 2, respectively. On the other hand, increasing the number of parameter (𝐹3 and 𝐹4), both conventional GA and PSO cannot get reasonable results especially for 𝐹3, while the proposed HPSO can perform much better results than conventional GA and PSO even for 𝐹3.

tab1
Table 1: Comparison results of test functions.

According to our experiments, we also found that our method is strongly affected by the number of subpopulations. Figure 5 shows the dependence of optimization accuracy on the number of subpopulations for 𝐹3. We change the number of subpopulations (𝑁) and perform 40 runs for each 𝑁. The averaged square error is shown in Figure 5. All of experiments use 5600 particles. We can see that the curve of average square error is getting down when the number of subpopulations (𝑁) is increased. The traditional PSO corresponds to 𝑁=1. This result indicates that our HPSO method can provide higher accuracy while using more subpopulation. However, if the number of subpopulations is too large, the accuracy will be decreased as shown in Figure 5 because of the limited particle number in a subpopulation.

561406.fig.005
Figure 5: The dependence of optimization accuracy on the number of subpopulations.
4.2. Medical Volume Data

In this part, we perform several registration experiments with medical volume phantom data to evaluate the performance of the proposed HPSO technique.

4.2.1. Simulated Data

We firstly use some simulated data to test our proposed method. Figure 6 shows some slices of our test data for 2 experiments which have different transformation parameters. The size of each volume is resized to 128×128×15, and the spacing of each volume is 2.59×2.59×8.0 (mm).

fig6
Figure 6: (a) MR volume (fixed image data); (b) CT volume (moving image data) of experiment 1 ((b) to (a)); (c) CT volume (moving image data) of experiment 2 ((c) to (a)).

These simulated data is made by specific parameters which we give it, so that we can easily estimate the result of registration accuracy by using these answers (actual results). Tables 2 and 3 show the comparative results of our HPSO and the previous method (GA and PSO) that diff is defined as difference between experimental results and ground truth.

tab2
Table 2: Comparison result of registration parameters ((b) to (a)).
tab3
Table 3: Comparison result of registration parameters ((c) to (a)).

According to the experimental results which are shown previously, we realize our proposed method can provide more exact result of registration parameters than conventional method. This is the first evidence to prove that our HPSO is a better optimization approach.

4.2.2. Vanderbilt Database

The real medical data (Vanderbilt database [11]) is also used to evaluate the performance of our proposed HPSO (see Figure 7). This database gives both multimodal brain volumes and their marker-based golden standard transforms. These rigid transforms are determined by marker-based prospective registration techniques and represented by eight couples of 3-D points (landmarks) on both of two medical volumes. Such a golden standard transform can be considered as a ground truth (correct one). The size of each volume is resized to 256×256×29 and the spacing of each volume is 1.25×1.25×4.0 (mm).

fig7
Figure 7: (a) A slice of CT image of Vanderbilt database; (b) a slice of MR image of Vanderbilt database.

Here we try gradient decent method, GA, PSO, and HPSO on CT to MR registration problem. A total of 3 runs with different random initial values were performed. Eight landmarks in the volume of Vanderbilt database are used for quantitative evaluation. Equation (12) is a measure used to compare the registration accuracy between experiment-obtained transform (𝑔) and the golden standard transform (𝑔Golden). 𝐱𝑖(𝑖=1,2,8) is coordinate of the landmark 𝑖: 1𝑒=88𝑖=1𝑔𝐱𝑖𝑔Golden𝐱𝑖.(13) The averaged accuracy of registration results for each method is shown in Table 4. It can be seen that our proposed HPSO performs much better results than conventional gradient decent method, GA and PSO.

tab4
Table 4: Comparison of registration accuracy (mm).
4.3. Parallel Implementation

The efficiency of registration is also an important issue as well as registration accuracy [12, 13]. Although we have shown the accuracy improvement in the previous sections, almost global optimization has a big disadvantage: a huge computation cost. We applied our proposed HPSO to a computer which has 4-core CPU and implement it with parallel technique to overcome this problem. For registering two images with size of 256×256×29, the average computation cost is reduced to 1893.637 sec from 3661.747 sec by our implementation. These results indicated that the computation cost can be significantly reduced by parallel implementation.

5. Conclusion and Future Works

This paper introduces a new global optimization approach named hybrid particle swarm optimization (HPSO) which incorporates two concepts: subpopulation and crossover of genetic algorithms into the conventional PSO. We performed both functional evaluation and 3-D rigid medical volume registration to estimate the proposed method. First, 4 functions have been applied to test the ability of our method for finding the global resolution and avoiding the local resolutions. Then, in case of 3-D rigid medical volume registration, we applied our HPSO to both simulated data and real medical data (Vanderbilt database) to discover the capability of this method for real medical volume registration. In order to compare conventional methods such as GA and PSO are also used for each experiment. All experimental results prove that the proposed HPSO performs much better results than conventional GA and PSO. Therefore, we can summarize that our HPSO is an advanced optimization method. The large computation cost can be significantly reduced by parallel implementation.

Furthermore, this work can be extended to 3D nonrigid registration in order to cover more computer-assisted surgery application such as real time liver tumor resection.

Acknowledgments

This work was supported in part by the Grant-in-Aid for Scientific Research from the Japanese Ministry for Education, Science, Culture and Sports under Grants no. 2430076 and no. 24103710 and in part by the R-GIRO Research fund from Ritsumeikan University.

References

  1. S. Morikawa, T. Inubushi, Y. Kurumi et al., “MR-guided microwave thermocoagulation therapy of liver tumors: Initial clinical experiences using a 0.5 T open MR system,” Journal of Magnetic Resonance Imaging, vol. 16, no. 5, pp. 576–583, 2002. View at Publisher · View at Google Scholar · View at Scopus
  2. R. Xu, Y. W. Chen, Y. T. Song, S. Morikawa, and Y. Kurumi, “Application of non-rigid medical image registration on open-MR based liver cancer surgery,” International Journal of Computer Assisted Radiology and Surgery, vol. 1, no. 7, pp. 286–288, 2006. View at Publisher · View at Google Scholar · View at Scopus
  3. L. G. Brown, “A survey of image registration techniques,” ACM Computing Surveys, vol. 24, no. 4, pp. 325–376, 1992. View at Publisher · View at Google Scholar · View at Scopus
  4. J. P. W. Pluim, J. B. A. Maintz, and M. A. Viergever, “Mutual-information-based registration of medical images: a survey,” IEEE Transactions on Medical Imaging, vol. 22, no. 8, pp. 986–1004, 2003. View at Publisher · View at Google Scholar · View at Scopus
  5. J. M. Rouet, J. J. Jacq, and C. Roux, “Genetic algorithms for a robust 3-D MR-CT registration,” IEEE Transactions on Information Technology in Biomedicine, vol. 4, no. 2, pp. 126–136, 2000. View at Publisher · View at Google Scholar · View at Scopus
  6. J. Kennedy and R. C. Eberhart, “Particle swarm optimization,” in Proceedings of the IEEE International Conference on Neural Networks, vol. 4, pp. 1942–1948, December 1995. View at Scopus
  7. Q. Bai, “Analysis of particle swarm optimization algorithm,” Computer and Information Science, vol. 3, no. 1, 2010.
  8. M. P. Wachowiak, R. Smolíková, Y. Zheng, J. M. Zurada, and A. S. Elmaghraby, “An approach to multimodal biomedical image registration utilizing particle swarm optimization,” IEEE Transactions on Evolutionary Computation, vol. 8, no. 3, pp. 289–301, 2004. View at Publisher · View at Google Scholar · View at Scopus
  9. Y.-W. Chen, C.-L. Lin, and A. Mimori, “Multimodal medical image registration using particle swarm optimization,” in Proceedings of the 8th International Conference on Intelligent Systems Design and Applications (ISDA '08), vol. 3, pp. 127–131, November 2008. View at Publisher · View at Google Scholar · View at Scopus
  10. D. Whitley, S. Rana, and K. Mathias, “Evaluating evolutionary algorithms,” Artificial Intelligence, vol. 85, pp. 245–276, 1996. View at Publisher · View at Google Scholar
  11. J. West, J. M. Fitzpatrick, M. Y. Wang et al., “Comparison and evaluation of retrospective intermodality brain image registration techniques,” Journal of Computer Assisted Tomography, vol. 21, no. 4, pp. 554–566, 1997. View at Publisher · View at Google Scholar · View at Scopus
  12. R. Shams, P. Sadeghi, R. A. Kennedy, and R. I. Hartley, “A survey of medical image registration on multicore and the GPU,” IEEE Signal Processing Magazine, vol. 27, no. 2, pp. 50–60, 2010. View at Publisher · View at Google Scholar · View at Scopus
  13. R. Shams, P. Sadeghi, R. Kennedy, and R. Hartley, “Parallel computation of mutual information on the GPU with application to real-time registration of 3D medical images,” Computer Methods and Programs in Biomedicine, vol. 99, no. 2, pp. 133–146, 2010. View at Publisher · View at Google Scholar · View at Scopus