Author name used in this publication: Michael FulhamAuthor name used in this publication: Dagan FengRefereed conference paper2008-2009 > Academic research: refereed > Refereed conference paperVersion of RecordPublishe

Eberl, S

Feng, DD

Fulham, MJ

Wen, L

Xia, Y

Author name used in this publication: Michael FulhamAuthor name used in this publication: Dagan FengRefereed conference paper2008-2009 > Academic research: refereed > Refereed conference paperVersion of RecordPublishedVoR allowe

PolyU Institutional Repository

Segmentation of Dual Modality Brain PET/CT Images Using the MAP-MRF Model Yong Xia*1, 3, Lingfeng Wen1, 2, 3, Stefan Eberl1, 2, Michael Fulham1, 2, 4, and Dagan Feng1, 31 Biomedical and Multimedia Information Technology (BMIT) Research Group School of Information Technologies, The University of Sydney, Australia 2 Department of PET and Nuclear Medicine, Royal Prince Alfred Hospital, Sydney, Australia 3 Department of Electronic & Information Engineering, Hong Kong Polytechnic University, Hong Kong 4 Faculty of Medicine, The University of Sydney, Australia y.xia@usyd.edu.au  Abstract—Dual modality PET/CT has now essentially replaced PET in clinical practice and provided an opportunity to improve image segmentation through the high resolution, lower noise CT data. Thus far most research efforts have concentrated on segmentation of PET-only data. In this work we propose a systematic solution for the automated segmentation of brain PET/CT images into gray, white matter and CSF regions with the MAP-MRF model. Our approach takes advantage of the full information available from the combined scan. A PET/CT image pair and its segmentation result are modelled as a random field triplet, and segmentation is eventually achieved by solving a maximum a posteriori (MAP) problem using the expectation-maximization (EM) algorithm with simulated annealing. We compared the novel algorithm to two widely used PET-only based segmentation methods in the SPM5 toolbox and the VBM toolbox for simulation and patient data. Our results suggest that using the proposed approach substantially improves the accuracy of the delineation of brain structures. I. INTRODUCTION Positron emission tomography (PET) can detect subtle functional changes at early stages of a disease process and this gives PET a distinct advantage over anatomical imaging techniques in the evaluation of neurodegenerative disorders [1]. However, accurate anatomic localization of functional abnormalities seen on PET scans is challenging [2] since PET data intrinsically have low spatial resolution, high noise and limited anatomical information. Recently PET has been combined with CT and the dual modality PET/CT scanner provides accurately aligned anatomic (CT) and functional (PET) images in a single scanning session. Although PET/CT has now essentially replaced PET in clinical practice, most functional brain image segmentation algorithms in the literature were designed for PET-only images [3, 4]. Among these algorithms, those in the SPM5 toolbox (Wellcome Department of Cognitive Neurology) [5] and the VBM Toolbox (University of Jena, Department of Psychiatry) [6] are two of the most widely used. The SPM algorithm is based on cluster analysis with a modified mixture model and prior information about the likelihood of each voxel belonging to each of three major brain structures. The VBM algorithm applies the Hidden Markov Random Field (HMRF) model and aims to minimize the noise level in the results. However, due to the lack of detailed anatomy, these algorithms rely heavily on the anatomical prior. The performance of these algorithms is greatly diminished when marked anatomical abnormalities seen with severe neurodegenerative disorders. We have previously developed a semi-automated technique to segment brain PET/CT images, in which the CSF region is delineated from CT data before gray matter and white matter is separated using both PET and CT data [7]. Our findings suggested that the high resolution, lower noise CT structural information improved segmentation performance. In this paper, we advance our previous work by proposing a systematic solution for the automated segmentation of brain PET/CT images by using a MAP-MRF model [8], which combines the Markov random field (MRF) theory with statistical decision and the maximum a posteriori (MAP) estimation [9]. We also compare the proposed approach to the SPM and VBM algorithms when applied to simulation and patient data. II. METHOD A. Simulation and Patient PET-CT Data The simulation was based on the Zubal anatomical phantom [12]. The PET sinogram data were simulated for an ECAT Exact HR+ scanner by the PET-SORTEO simulator [13] assuming FDG activity of 23.0 MBq in gray matter and 8.5 MBq in white matter. A 3D filtered back projection method was applied to derive reconstructed PET data with corrections for random, dead-time, attenuation and scatter. The CT data were simulated by assuming (1) initial values for the background, gray matter, white matter, CSF, skull, and non-brain tissues of  -986.6, 61.7, 28.3, 11.2, 1027.0, and 23.1 Hounsfield Units (HU), respectively based on clinical data; (2) the point spread function (PSF) is a Gaussian function (FWHM = [4.4, 4.4, 4.2] mm); and (3) the noise in these regions was assumed to be Gaussian noise with standard deviations of 69.6, 7.7, 5.6, 8.9 364.2, and 9.1 HU, respectively, again based on clinical data [7]. The simulated data have dimensions of 128256256 ××  and voxel sizes of 4.11.11.1 ××  mm. 978-1-4244-2295-1/08/$25.00 © 2008 IEEE MMSP 2008107Authorized licensed use limited to: Hong Kong Polytechnic University. Downloaded on June 24,2010 at 07:28:01 UTC from IEEE Xplore.  Restrictions apply. Clinical data were selected from the clinical studies from the Department of PET and Nuclear Medicine at Royal Prince Alfred Hospital and were acquired on a Biograph LSO Duo PET/CT scanner. Each clinical data set has a dimension of 47256256 ××  and a voxel size of 4.3977.0977.0 ××  mm. B. Background Removal Before segmentation, background, skull and non-brain tissues were removed from each brain PET/CT image through the following steps. i) The brain-mask provided by the SPM5 toolbox [14] was spatially normalized to the PET data and then dilated to ensure all brain tissues were covered by it. ii) The resultant brain-mask was applied to the CT data to remove background, non-brain tissues and most of the skull. iii) The remnant CT data consisted of brain and some skull. Since the HU values for bones are usually higher than 200, the boundary of the brain in CT data was refined by removing all voxels whose HU values were outside the range [-10, 70]. iv) Morphological operators were finally used to close holes inside the brain regions. After background removal, only the voxels belonging to brain regions were subsequently processed by the segmentation algorithm. C. MAP-MRF Model Each slice of a brain PET/CT image pair, together with its segmentation result, which is a set of labels and is called a configuration in the terminology of MRF theory, was considered as a random field triplet LGG CTPET ,,  defined on a rectangular lattice ( ){ }HjWijiS ≤≤≤≤= 1,1:, . An observed PET image { }Ssgg PETsPET ∈= : , a CT image { }Ssgg CTsCT ∈= :  and a configuration { }Ssll s ∈= :  are instances of random field PETG , CTG  and L , respectively. According to the MAP theory [8, 9], the segmentation process is equivalent to the search for the optimal configuration *l , which maximizes the posterior probability on condition of the observed PET-CT image. With the Bayes rule, it can be formalized as ( ) ( )( )CTPETCTPETlgg ggplplggplLCTPET,,max|* , Ω∈= ,                (1) where LΩ  is the configuration space, ( )lggp CTPET ,  is the likelihood probability of the observed PET-CT image on condition of the configuration l , ( )lp  is the prior probability of the configuration l , and ( )CTPET ggp ,  is the probability of the observed image, which does not vary with respect to any configurations. In this paper, the probability ( )xXp =  is denoted as ( )xp  for simplicity. Assuming the PET random field PETG  is independent of the CT random field CTG  with respect to a given configuration l , the Eq. (1) can be simplified as ( ) ( ) ( )lplgplgpl CTPETlgg LCTPETΩ∈= max|* , .                (2) Two assumptions were made to estimate the probabilities in Eq. (2). One is that the voxel values of each targeted region in either image follow an independent and identical Gaussian distribution. Consequently, the likelihood probability can be calculated as ( ) ( )( )∏ ⎟⎟⎠⎞⎜⎜⎝⎛−−=+++++s ccscglgp 222exp21σπμσπ,              (3) where “+” represents either PET or CT, and ( )++ cc σμ ,  are the distribution parameters. The other assumption is that the configuration field L  is a second order MRF [9] and its joint distribution function is the following Gibbs function [15] ( ) ( )⎥⎥⎦⎤⎢⎢⎣⎡−= ∑∑∈s ttssllVTZlPη,1exp1 ,                       (4) where Z  is a normalizing constant, T  is a constant analogous to temperature, sη  is the second order neighbourhood of s , and ( )ts llV ,  is the potential of the pair-wised clique { }ts ll , , which is defined as ( )⎩⎨⎧=−≠=tststs llllllV,1,1, .                          (5) With these assumptions, the MAP-MRF model-based segmentation of PET-CT images can be interpreted as minimizing the following energy function ( )[ ]∑ ++= Ω∈ s LsCTsPETslgg EEEl LCTPET αmin|* , ,             (6) where ( )( ) ( ) ( )πσσπμ 2ln21ln2 22++−=++++=+cccsclsgEs              (7) is the Gaussian energy, ( )∑∈=sttsLs llVEη,                                 (8) is the Gibbs energy, and α  is a weighting factor. D. Segmentation The segmentation problem as formulated in Eq. (7) was solved by the expectation-maximization (EM) algorithm with simulated annealing [10, 11], as outlined below. 1. Initialization: Half of the labels were randomly initialized, and other labels were initialized according to the prior of gray matter, white matter, and CSF in the SPM5 toolbox. 2. E-step: Use the maximum likelihood estimation [10] to calculate Gaussian distribution parameters. 3. M-step: Use the Metropolis sampler [8] with the simulated annealing scheme [16] to refine the configuration. The labels were updated in a raster scan order. For each voxel s , a new configuration 'l  was created by randomly choosing a new label for that voxel. The energy of these two configurations can be evaluated by comparing only the following terms clLsclCTsclPETsclsssssEEEE====+⎟⎠⎞⎜⎝⎛ += α               (9) '''''clLsclCTsclPETsclsssssEEEE====+⎟⎠⎞⎜⎝⎛ += α .         (10) 108Authorized licensed use limited to: Hong Kong Polytechnic University. Downloaded on June 24,2010 at 07:28:01 UTC from IEEE Xplore.  Restrictions apply.                  (a)                                (b)                                (c)                                (d)                                 (e)                                (f)                               (g) Fig. 1. The 31st slice of a clinical study and its segmentation results. (a) PET image, (b) CT image, (c) CT image (after background removal), (d) Initial configuration, (e) Configuration after 10 iterations, (f) Configuration after 20 iterations, (g) Final configuration (after 32 iterations).                      (a)                                      (b)                                     (c)                                       (d)                                       (e)                                     (f) Fig. 2. The 74th slice of the PET-CT simulation data and its segmentation results. (a) PET image, (b) CT image (after background removal), (c) Atlas, (d) Result of the SPM algorithm, (e) Result of the VBM algorithm, (f) Result of the proposed algorithm.  If ss EE <' , accept the new configuration 'l . Otherwise, the new configuration 'l  will be accepted with the probability ( ) ( ) ⎟⎟⎠⎞⎜⎜⎝⎛ −−=nTEEnp sss'exp ,                     (11) where  ( ) ( )10 += nTnT                              (12) is the annealing temperature sequence, 10 =T  is the initial temperature, and n  is the iteration number. 4. Repeat above two steps until i) the number of iterations exceeds 50; or ii) less than 15 voxels changed their labels in the latest iteration. The weighting factor α  presents the trade-off between the fidelity of the observed image and the smoothness prior of the configuration. It is expected to be relatively large in the early stage of the optimization so that the voxel values can play a dominant role, and considerably smaller in the later stage to emphasise the spatial constraints. To this end, the following attenuation function is applied to the variable α . ( ) 10 9.0 ααα +×= nn ,                                  (13) where 10 αα +  is the initial weight, and the final weight approaches 1α . In this study, we empirically set 300 =α  and 11 =α . III. RESULTS The 31st slice of a clinical study and its segmentation results are shown in Fig. 1, where the voxel values of PET and CT data are linearly mapped to [0, 255] for display. Comparing Fig. 1(b) and (c), it is recognized that the contrast of three major brain structures has been markedly improved in background-free CT data, as a result of the dramatic decrease in dynamic range of the CT data after background removal.  In Fig. 1(e)-(g), it is demonstrated that the proposed approach sacrifices preservation of certain detail for the ability to achieve a better noise resistance during the search for an optimal configuration.  Next, the proposed approach was compared to the SPM and the VBM algorithms for the PET/CT simulation data. The performance of these algorithms was assessed with the overlapping area of the segmentation results and the atlas, which was quantitatively measured by the Dice similarity coefficient (DSC) [17] AtlasSegAtlasSegAAAADSC+∩= 2 ,                          (14) where SegA  is the segmented area and AtlasA  is the atlas area. Generally, a better segmentation is indicated by a higher DSC value. The segmentation results of the 74th slice are shown in Fig. 2. The DSC values calculated on the results of nine slices are depicted in Fig. 3. It reveals that, though not performing best in the gray matter delineation, the proposed approach substantially outperforms the other methods in accurate delineation of white matter and CSF.                          (a)                                     (b)                                     (c) Fig. 3. Performance of three segmentation approaches on 9 slices of the PET-CT simulation data. (1) Delineation of gray matter, (2) Delineation of white matter, (3) Delineation of CSF. 109Authorized licensed use limited to: Hong Kong Polytechnic University. Downloaded on June 24,2010 at 07:28:01 UTC from IEEE Xplore.  Restrictions apply.                          (a)                                               (b)                                               (c)                                               (d)                                             (e) Fig. 4. The 29th slice (top row) and the 27th slice (bottom row) from two clinical studies and their segmentation results. (a) PET image, (b) CT image (after background removal), (c) Result of the SPM algorithm, (d) Result of the VBM algorithm, (e) Result of the proposed algorithm.  Finally, the comparative experiments were carried out on several clinical studies. The results acquired on the 29th and 27th slices of two clinical studies are shown in Fig. 4. It illustrates again that our approach provides more accurate brain PET/CT image segmentation.  IV. CONCLUSION We propose an automated segmentation algorithm for brain PET/CT images based on the MAP-MRF model. By combining PET and CT information, our image segmentation approach provides a more accurate delineation of gray, white matter and CSF regions than PET-only based methods. Our technique requires minimal user interaction and supports the great potential of using the information from multiple modalities for robust and effective segmentation of medical images. Our future work is aimed at extending the technique to the segmentation of 3D volume PET/CT data and more efficient methods to integrate PET and CT information.  ACKNOWLEDGMENT This work was supported in part by ARC and PolyU grants. REFERENCES [1] J.R. Petrella, R.E. Coleman, and P. M. Doraiswamy, “Neuroimaging and early diagnosis of Alzheimer Disease: A look to the Future,” Radiology, vol. 226, no. 2, pp. 315-336, Feb. 2003. [2] D.W. Townsend, J.P.J. Carney, J.T. Yap, and N.C. Hall, “PET/CT today and tomorrow”, J. Nucl. Med., vol. 45, pp. 4S-14S, 2004. [3] J. Mykkänen, M. Juhola, and U. Ruotsalainen, “Extracting VOIs from brain PET images,” International Journal of Medical Informatics, vol. 58-59, pp. 51-57, Sep. 2000.  [4] W. Zhu, T. Jiang, “Automation segmentation of PET image for brain tumors,” in Conf. Rec. 2003 IEEE Conf. Nuclear Science Symp., vol. 4, pp. 2627–2629. [5] J. Ashburner, and K.J. Friston, “Multimodal Image Coregistration and Partitioning - a Unified Framework,” NeuroImage, vol. 6, pp. 209-217, 1997. [6] J. Ashburner, and K.J. Friston, “Unified segmentation,” Neuroimage, vol. 26, no. 3, pp. 839-851, 2005. [7] Y. Xia, L. Wen, S. Eberl, M. Fulham, and D. Feng, "Segmentation of brain structures using PET-CT images," in Proc. of the 5th International Conference on Information Technology and Applications in Biomedicine (ITAB 2008), Shenzhen, China, May 30-31, 2008. [8] S. Geman, and D. Geman, “Stochastic relaxation, Gibbs distributions and the Bayesian restoration of images,” IEEE Trans. Pattern Anal. Machine Intell., vol. MPAI-6, no. 6, pp. 721–741, Nov. 1984.  [9] S.Z. Li, “Modeling image analysis problems using Markov random fields,” in C. R. Rao and D. N. Shanbhag (ed), Stochastic Processes: Modeling and simulation, vol. 20 of Handbook of Statistics. Elsevier Science, 2000, pp. 1-43.  [10] H. Deng, and D.A. Clausi, “Unsupervised image segmentation using a simple MRF model with a new implementation scheme,” Pattern Recognit., vol. 37, no. 12, pp. 2323-2335, 2004. [11] Y. Xia, D. Feng, and R. Zhao, “Adaptive Segmentation of Textured Images by Using the Coupled Markov Random Field Model,” IEEE Trans. Image Processing, vol. 15, pp. 3559–3566, Nov. 2006. [12] I. G. Zubal, C. R. Harrell, E. O. Smith, Z. Rattner, G. Gindi, and P.B. Hoffer, “Computerized 3-dimensional segmented human anatomy,” Med. Phys., vol. 21, pp. 299–302, Feb. 1994. [13] A. Reilhac, G. Batan, C. Michel, C. Grova, J. Tohka, N. costes, and A.C. Evans, “PET-SORTEO: Validation and Development of database of simulated PET volumes,” IEEE Trans. on Nucl. Sci., vol. 52, no. 5, pp. 46-52, 2005. [14] R. S. J. Frackowiak, K. J. Friston, C. D. Frith, et al., Human Brain Function. Amsterdam; Boston: Elsevier Academic Press, 2004. [15] J.M. Hammersley, and P. Clifford, “Markov field on finite graphs and lattices”, unpublished, 1971.  [16] S. Kirkpatrick, C. D. Gelatt Jr, and M. P. Vecchi, “Optimization by simulated annealing,” Science, vol. 220, no. 4598, pp. 671-680, 1983. [17] A.P. Zijdenbos, B.M. Dawant, R.A. Margolin, and A.C. Palmer, “Morphometric analysis of white matter lesions in MR images: Methods and validation,” IEEE Trans. on Med. Imag., vol. 13, no. 4, pp. 716–724, Dec. 1994.  110Authorized licensed use limited to: Hong Kong Polytechnic University. Downloaded on June 24,2010 at 07:28:01 UTC from IEEE Xplore.  Restrictions apply. 

Segmentation of dual modality brain PET/CT images using the MAP-MRF model

Abstract

Similar works

Full text

Available Versions

PolyU Institutional Repository