Abstract:
A method for segmenting an anatomical image, including: receiving a patient anatomical image; receiving a baseline anatomical image having pre-segmented labels, wherein the pre-segmented labels identify regions of interest in the baseline anatomical image; aligning the patient anatomical image with the baseline anatomical image to produce a transformation that when applied to the pre-segmented labels roughly identifies regions of interest in the patient anatomical image that correspond to the regions of interest in the baseline anatomical image; and updating the pre-segmented labels, which have been deformed by application of the transformation, with a new transformation that minimizes the likelihood of intensity distributions within the regions of interest of the patient anatomical image to produce a gradient image that better identifies the regions of interest of the patient anatomical image.

Description:
CROSS-REFERENCE TO RELATED APPLICATION 
       [0001]    This application claims the benefit of U.S. Provisional Application No. 61/032,237, filed Feb. 28, 2008, the disclosure of which is incorporated by reference herein in its entirety. 
     
    
     BACKGROUND OF THE INVENTION 
       [0002]    1. Technical Field 
         [0003]    The present invention relates to multi-label segmentation and segmenting different organs of the abdomen. 
         [0004]    2. Discussion of the Related Art 
         [0005]    Image segmentation is the process of partitioning an image into different regions. A goal of image segmentation is to obtain a higher-level description of image content. For instance, in medical imaging, the segmentation of anatomical structures is a key element for computer-aided diagnosis and image-guided therapies. 
       SUMMARY OF THE INVENTION 
       [0006]    In an exemplary embodiment of the present invention, a method for segmenting an anatomical image, comprises: receiving a patient anatomical image: receiving a baseline anatomical image having pre-segmented labels, wherein the pre-segmented labels identify regions of interest in the baseline anatomical image; aligning the patient anatomical image with the baseline anatomical image to produce a transformation that when applied to the pre-segmented labels roughly identifies regions of interest in the patient anatomical image that correspond to the regions of interest in the baseline anatomical image; and updating the pre-segmented labels, which have been deformed by application of the transformation, with a new transformation that minimizes the likelihood of intensity distributions within the regions of interest of the patient anatomical image to produce a gradient image that better identifies the regions of interest of the patient anatomical image. 
         [0007]    The method further comprises computing the new transformation, wherein computing the new transformation comprises: computing a gradient for all the regions of interest of the patient anatomical image; regularizing the gradient; and generating the new transformation by using the regularized gradient. 
         [0008]    The new transformation is applied to the deformed pre-segmented labels by computing a composition of the deformed pre-segmented labels and the new transformation. 
         [0009]    Computing the gradient for all the regions of interest of the patient anatomical image comprises: (1) for a region of interest of the patient anatomical image, computing a temporary image for the region of interest; computing an intensity distribution for the region of interest; and computing a gradient for the region of interest; (2) updating the gradient image with the gradient for the region of the interest; and repeating (1) and (2) until the gradient image has been updated with a gradient for all the regions of interest of the patient anatomical image. 
         [0010]    The pre-segmented labels are repeatedly updated with new transformations until all the regions of interest of the patient anatomical image are better identified. 
         [0011]    The patient anatomical image comprises an abdomen. 
         [0012]    The patient anatomical image is a computed tomography (CT) image. 
         [0013]    In an exemplary embodiment of the present invention, a system for segmenting an anatomical image, comprises: a memory device for storing a program: a processor in communication with the memory device, the processor operative with the program to: receive a patient anatomical image; receive a baseline anatomical image having pre-segmented labels, wherein the pre-segmented labels identify regions of interest in the baseline anatomical image; align the patient anatomical image with the baseline anatomical image to produce a transformation that when applied to the pre-segmented labels roughly identifies regions of interest in the patient anatomical image that correspond to the regions of interest in the baseline anatomical image; and update the pre-segmented labels, which have been deformed by application of the transformation, with a new transformation that minimizes the likelihood of intensity distributions within the regions of interest of the patient anatomical image to produce a gradient image that better identifies the regions of interest of the patient anatomical image. 
         [0014]    The processor is further operative with the program to compute the new transformation, wherein when computing the new transformation the processor is further operative with the program to: compute a gradient for all the regions of interest of the patient anatomical image; regularize the gradient; and generate the new transformation by using the regularized gradient. 
         [0015]    The new transformation is applied to the deformed pre-segmented labels by computing a composition of the deformed pre-segmented labels and the new transformation. 
         [0016]    When computing the gradient for all the regions of interest of the patient anatomical image the processor is further operative with the program to: (1) for a region of interest of the patient anatomical image, compute a temporary image for the region of interest; compute an intensity distribution for the region of interest; and compute a gradient for the region of interest; (2) update the gradient image with the gradient for the region of the interest; and repeat (1) and (2) until the gradient image has been updated with a gradient for all the regions of interest of the patient anatomical image. 
         [0017]    The pre-segmented labels are repeatedly updated with new transformations until all the regions of interest of the patient anatomical image are better identified. 
         [0018]    The patient anatomical image comprises an abdomen. 
         [0019]    The patient anatomical image is a CT image. 
         [0020]    In an exemplary embodiment of the present invention, a computer readable medium tangibly embodying a program of instructions executable by a processor to perform method steps for segmenting an anatomical image is provided, the method steps comprising: receiving a patient anatomical image; receiving a baseline anatomical image having pre-segmented labels, wherein the pre-segmented labels identify regions of interest in the baseline anatomical image; aligning the patient anatomical image with the baseline anatomical image to produce a transformation that when applied to the pre-segmented labels roughly identifies regions of interest in the patient anatomical image that correspond to the regions of interest in the baseline anatomical image; and updating the pre-segmented labels, which have been deformed by application of the transformation, with a new transformation that minimizes the likelihood of intensity distributions within the regions of interest of the patient anatomical image to produce a gradient image that better identifies the regions of interest of the patient anatomical image. 
         [0021]    The method steps further comprise computing the new transformation, wherein computing the new transformation comprises: computing a gradient for all the regions of interest of the patient anatomical image; regularizing the gradient; and generating the new transformation by using the regularized gradient. 
         [0022]    The new transformation is applied to the deformed pre-segmented labels by computing a composition of the deformed pre-segmented labels and the new transformation. 
         [0023]    Computing the gradient for all the regions of interest of the patient anatomical image comprises: (1) for a region of interest of the patient anatomical image, computing a temporary image for the region of interest; computing an intensity distribution for the region of interest; and computing a gradient for the region of interest; (2) updating the gradient image with the gradient for the region of the interest; and repeating (1) and (2) until the gradient image has been updated with a gradient for all the regions of interest of the patient anatomical image. 
         [0024]    The pre-segmented labels are repeatedly updated with new transformations until all the regions of interest of the patient anatomical image are better identified. 
         [0025]    The patient anatomical image comprises an abdomen. 
         [0026]    The patient anatomical image is a CT image. 
         [0027]    The foregoing features are of representative embodiments and are presented to assist in understanding the invention. It should be understood that they are not intended to be considered limitations on the invention as defined by the claims, or limitations on equivalents to the claims. Therefore, this summary of features should not be considered dispositive in determining equivalents. Additional features of the invention will become apparent in the following description, from the drawings and from the claims. 
     
    
     
       BRIEF DESCRIPTION OF THE DRAWINGS 
         [0028]      FIGS. 1A-C  are images that illustrate multi-label segmentation according to an exemplary embodiment of the present invention; 
           [0029]      FIGS. 2A  and B are flowcharts that illustrate a method for multi-label segmentation according to an exemplary embodiment of the present invention; and 
           [0030]      FIG. 3  is a block diagram of a system in which exemplary embodiments of the present invention may be implemented. 
       
    
    
     DETAILED DESCRIPTION OF EXEMPLARY EMBODIMENTS 
       [0031]    A hierarchical multi-label segmentation method based on non-rigid registration techniques to segment an arbitrary number of regions, according to an exemplary embodiment of the present invention, will hereinafter be described. In an exemplary embodiment of the method, first align an image I S , with pre-segmented labels I T     N   , to the image to be segmented I. Then, deform the pre-segmented labels I T     N    and use them as a rough initialization to a multi-label segmentation technique, according to an exemplary embodiment of the present invention, where the deformed pre-segmented labels I T     N   , are non-rigidly aligned to the image I by maximizing the likelihood of intensity distributions within different regions of interest. The intensity models and the corresponding posteriori distributions are estimated and updated throughout the alignment. The method according to an exemplary embodiment of the present invention allows a spatial relation between different regions of interest to be kept by finding local variations of shapes through one deformation field. An example of the method according to an exemplary embodiment of the present invention applied to segment eight regions of computed tomography (CT) images of the abdomen, is further described hereinafter. 
         [0032]    A description of the statistical formulation of region-based segmentation will now be provided. 
         [0033]    Let Ω ε R d  be open and bounded, and I:Ω→R be the image to be segmented. Assume that Ω is a partition composed of N independent disjoint regions Ω i . This gives the simplified expression: 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       p 
                        
                       
                         ( 
                         
                           I 
                           | 
                           
                             P 
                              
                             
                               ( 
                               Ω 
                               ) 
                             
                           
                         
                         ) 
                       
                     
                     = 
                     
                       
                         p 
                          
                         
                           ( 
                           
                             I 
                             | 
                             
                               { 
                               
                                 
                                   Ω 
                                   1 
                                 
                                 , 
                                 … 
                                  
                                 
                                     
                                 
                                 , 
                                 
                                   Ω 
                                   N 
                                 
                               
                               } 
                             
                           
                           ) 
                         
                       
                       = 
                       
                         
                           ∏ 
                           
                             i 
                             = 
                             1 
                           
                           N 
                         
                          
                         
                             
                         
                          
                         
                           p 
                            
                           
                             ( 
                             
                               I 
                               | 
                               
                                 Ω 
                                 i 
                               
                             
                             ) 
                           
                         
                       
                     
                   
                   , 
                 
               
               
                 
                   ( 
                   1 
                   ) 
                 
               
             
           
         
       
     
         [0000]    where p(I|Ω i ) denotes the probability of the image I where Ω i  is the region of interest. Assume that values of I at different locations of the same region can be modeled as an independent and identically distributed realization of the same random process. Define p i (I(x)) as the probability density function of a random variable modeling intensity values I(x) in Ω i . Given this model, the optimal partition can be obtained using a maximum likelihood principle, and minimizing the following energy proposed in [Zhu, S. C., Yuille, A. L.: Region competition: Unifying snakes, region growing, and bayes/MDL for multiband image segmentation. IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 18(9), 1996, pp. 884-900], the disclosure of which is incorporated by reference herein in its entirety: 
         [0000]    
       
         
           
             
               
                 
                   
                     E 
                      
                     
                       ( 
                       
                         { 
                         
                           
                             Ω 
                             i 
                           
                           , 
                           … 
                            
                           
                               
                           
                           , 
                           
                             Ω 
                             N 
                           
                         
                         } 
                       
                       ) 
                     
                   
                   = 
                   
                     
                       ∑ 
                       
                         i 
                         = 
                         1 
                       
                       N 
                     
                      
                     
                       
                         ∫ 
                         
                           Ω 
                           i 
                         
                       
                        
                       
                         
                           - 
                           log 
                         
                          
                         
                             
                         
                          
                         
                           
                             p 
                             i 
                           
                            
                           
                             ( 
                             
                               I 
                                
                               
                                 ( 
                                 x 
                                 ) 
                               
                             
                             ) 
                           
                         
                          
                         
                             
                         
                          
                         
                           
                              
                             x 
                           
                           . 
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   2 
                   ) 
                 
               
             
           
         
       
     
         [0034]    In the context of contour evolution, this energy can be expressed as the following energy to minimize: 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       E 
                        
                       
                         ( 
                         
                           
                             Ω 
                             i 
                           
                           , 
                           
                             p 
                             i 
                           
                         
                         ) 
                       
                     
                     = 
                     
                       
                         ∑ 
                         
                           i 
                           = 
                           1 
                         
                         N 
                       
                        
                       
                         ( 
                         
                           
                             - 
                             
                               
                                 ∫ 
                                 
                                   Ω 
                                   i 
                                 
                               
                                
                               
                                 log 
                                  
                                 
                                     
                                 
                                  
                                 
                                   p 
                                   i 
                                 
                               
                             
                           
                           - 
                           
                             
                               v 
                               2 
                             
                              
                             
                               
                                 ∫ 
                                 
                                   Γ 
                                   i 
                                 
                               
                                
                               
                                   
                               
                                
                               
                                  
                                 s 
                               
                             
                           
                         
                         ) 
                       
                     
                   
                   , 
                 
               
               
                 
                   ( 
                   3 
                   ) 
                 
               
             
           
         
       
     
         [0000]    where Γ i  represent the contour of the region Ω i , and the parameter ν controls the length of the contours. In particular, this energy is expressed in the context of level sets with a function φ i  that represents the region Ω i  where φ i (x)&gt;0 if and only if x ε Ω i : 
         [0000]    
       
         
           
             
               
                 
                   
                     E 
                      
                     
                       ( 
                       
                         
                           ϕ 
                           i 
                         
                         , 
                         
                           p 
                           i 
                         
                       
                       ) 
                     
                   
                   = 
                   
                     
                       ∑ 
                       
                         i 
                         = 
                         1 
                       
                       N 
                     
                      
                     
                       ( 
                       
                         
                           - 
                           
                             
                               ∫ 
                               Ω 
                             
                              
                             
                               
                                 H 
                                  
                                 
                                   ( 
                                   
                                     ϕ 
                                     i 
                                   
                                   ) 
                                 
                               
                                
                               log 
                                
                               
                                   
                               
                                
                               
                                 p 
                                 i 
                               
                                
                               
                                  
                                 x 
                               
                             
                           
                         
                         + 
                         
                           
                             v 
                             2 
                           
                            
                           
                              
                             
                               vH 
                                
                               
                                 ( 
                                 
                                   ϕ 
                                   i 
                                 
                                 ) 
                               
                             
                              
                           
                            
                           
                              
                             x 
                           
                         
                       
                       ) 
                     
                   
                 
               
               
                 
                   ( 
                   4 
                   ) 
                 
               
             
           
         
       
     
         [0035]    This formulation does not respect implicitly the condition of disjoint regions, but the minimization of this energy ensures that a pixel is assigned to only one region according to the maximum likelihood principle. 
         [0036]    A description of the method for non-rigid registration according to an exemplary embodiment of the present invention will now be provided. 
         [0037]    In the following description, given two images I 1  and I 2 , the registration problem is formulated as finding a mapping φ:Ω→Ω that maximizes a similarity measure between the images: S(I 1 , I 2 ∘φ). First, maximize the local cross correlation between I and I S , S LCC (I,I T     N   ∘φ) and apply the mapping φ to I T     N   . Second, maximize the likelihood of intensity distributions within different regions of interest: the multi-label similarity measure S ML (I,I T     N   ∘φ). This similarity measure according to an exemplary embodiment of the present invention allows the segmentation of different regions of interest to be refined. 
         [0038]    To find the optimal high-dimensional transformation, a sequence of transformations (φ k )k=0, . . . ,+∞, is built by composition of small displacements as described in [Chefd&#39;hotel, C., Hermosillo, G., Faugeras, O.: Flows of diffeomorphisms for multimodal image registration. In: Proceedings of IEEE International Symposium on Biomedical Imaging. (2002), pp. 753-756], the disclosure of which is incorporated by reference herein in its entirety, 
         [0000]      φ k+1 =φ k ∘(φ id +αν k ), φ 0 =φ id ,   (5) 
         [0000]    where φ id  is the identity transformation and ν k  is a velocity vector field that follows the gradient of the cost function to be minimized. Here, ν k  is obtained by computing the variational gradient of the cost function of the Local Cross-Correlation (LCC) similarity measure, i.e., ∇S LCC (I,I S ∘φ) or the ML similarity measure ∇S ML (I,I T     N   ∘φ). 
         [0039]    The gradient ν k  is regularized using a fast recursive filtering technique. This approximates a Gaussian smoothing, as described, for example, in [Deriche, R.: Recursively implementing the Gaussian and its derivatives. In: Proceedings of the International Conference on Image Processing, Singapore (September 1992), pp. 263-267], that has proven very efficient in practice. Here, deriving the similarity measure energy according to a high-dimensional transformation results in a vector field ν. To guarantee a well-posed problem, this vector field has to be regularized. For this purpose, different techniques have been proposed. The approach proposed in [Christensen, G. E., Rabbit, R. D., Miller, M. I.: Deformable templates using large deformation kinematics. IEEE Transactions on Image Processing, vol. 5(10), 1996, pp. 1435-1447], the disclosure of which is incorporated by reference herein in its entirety, solves the registration problem using a partial differential equation and has the advantage of capturing large deformations. In the method according to an exemplary embodiment of the present invention, a Gaussian filtering is used that can be seen as a variant of the fluid-approach described in Christensen et al. 
         [0040]    The previous iterative scheme (Eq. 5) is repeated until convergence, and can be seen as the discretization (via Taylor expansion) of the transport equation in the Eulerian frame: 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       
                         ∂ 
                         
                           φ 
                           t 
                         
                       
                       
                         ∂ 
                         t 
                       
                     
                     = 
                     
                       
                         - 
                         D 
                       
                        
                       
                           
                       
                        
                       
                         
                           φ 
                           t 
                         
                         · 
                         v 
                       
                     
                   
                   , 
                   
                     
                       φ 
                       0 
                     
                     = 
                     
                       φ 
                       id 
                     
                   
                   , 
                 
               
               
                 
                   ( 
                   6 
                   ) 
                 
               
             
           
         
       
     
         [0000]    where Dφ t  stands for the Jacobian matrix of φ t . Here, large deformations are possible because the regularization is applied to the velocity rather than the deformation described in [Dupuis, P., Grenander, U., Miller, M.: Variational problems on flows of diffeomorphisms for image matching. Quarterly of Applied Mathematics LVI(3), (1998), pp. 587-600], which details the suitable regularity conditions on the velocity field to generate a diffeomorphism. 
         [0041]    The method according to an exemplary embodiment of the present invention is embedded in a coarse-to-fine strategy. This reduces the computational cost by working with less data at lower resolutions. This also allows large displacements to be recovered, and helps avoiding local minima. In the method according to an exemplary embodiment of the present invention, five-levels of multi-resolutions are used. 
         [0042]    To refine the segmentation, in accordance with an exemplary embodiment of the present invention, a multi-labeled template matching algorithm that recovers local deformations of the shape obtained in the previous section is provided. Consider the registration framework, an image I T     N    ε Ω composed of N disjoint regions is defined, each region with a different label. This image can be seen as the union of N images representing a different region: 
         [0000]    
       
         
           
             
               
                 
                   
                     I 
                     
                       T 
                       N 
                     
                   
                   ≅ 
                   
                     
                       ∑ 
                       
                         i 
                         = 
                         1 
                       
                       N 
                     
                      
                     
                       
                         ( 
                         Ii 
                         ) 
                       
                       . 
                     
                   
                 
               
               
                 
                   ( 
                   7 
                   ) 
                 
               
             
           
         
       
     
         [0043]    Formulate the problem as finding a transformation φ:Ω→Ω that minimizes the likelihood between the intensity distribution functions of different regions p i  according to I and I T     N   . Thus, the following energy is minimized: 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       S 
                       ML 
                     
                      
                     
                       ( 
                       
                         I 
                         , 
                         
                           
                             I 
                             
                               T 
                               N 
                             
                           
                           · 
                           φ 
                         
                       
                       ) 
                     
                   
                   = 
                   
                     - 
                     
                       
                         ∫ 
                         Ω 
                       
                        
                       
                         
                           ( 
                           
                             
                               ∑ 
                               
                                 i 
                                 = 
                                 1 
                               
                               N 
                             
                              
                             
                               
                                 ( 
                                 
                                   
                                     I 
                                     i 
                                   
                                   · 
                                   
                                       
                                   
                                    
                                   φ 
                                 
                                 ) 
                               
                                
                               log 
                                
                               
                                   
                               
                                
                               
                                 
                                   p 
                                   i 
                                 
                                  
                                 
                                   ( 
                                   
                                     I 
                                      
                                     
                                       ( 
                                       x 
                                       ) 
                                     
                                   
                                   ) 
                                 
                               
                             
                           
                           ) 
                         
                          
                         
                           
                              
                             x 
                           
                           . 
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   8 
                   ) 
                 
               
             
           
         
       
     
         [0044]    In this equation, I T     N   ∘φ is the warped multi-labeled template and ∘ the composition operator. Since an optimal transformation φ is wanted, the derivation of the energy leads to the following gradient descent: 
         [0000]    
       
         
           
             
               
                 
                   
                     ∇ 
                     
                       
                         S 
                         ML 
                       
                        
                       
                         ( 
                         
                           I 
                           , 
                           
                             
                               I 
                               
                                 T 
                                 N 
                               
                             
                             · 
                             φ 
                           
                         
                         ) 
                       
                     
                   
                   = 
                   
                     
                       ∑ 
                       
                         i 
                         = 
                         1 
                       
                       N 
                     
                      
                     
                       
                         ∇ 
                         
                           ( 
                           
                             
                               I 
                               i 
                             
                             · 
                             φ 
                           
                           ) 
                         
                       
                        
                       
                         
                           ( 
                           
                             log 
                              
                             
                                 
                             
                              
                             
                               
                                 p 
                                 i 
                               
                                
                               
                                 ( 
                                 
                                   I 
                                    
                                   
                                     ( 
                                     x 
                                     ) 
                                   
                                 
                                 ) 
                               
                             
                           
                           ) 
                         
                         . 
                       
                     
                   
                 
               
               
                 
                   ( 
                   9 
                   ) 
                 
               
             
           
         
       
     
         [0045]    The density probability function of different regions is as follows: 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       p 
                       i 
                     
                      
                     
                       ( 
                       j 
                       ) 
                     
                   
                   = 
                   
                     
                       1 
                       
                          
                         
                           Ω 
                           i 
                         
                          
                       
                     
                      
                     
                       
                         ∫ 
                         Ω 
                         
                             
                         
                       
                        
                       
                         
                           ( 
                           
                             
                               I 
                               i 
                             
                             · 
                             φ 
                           
                           ) 
                         
                          
                         
                           ∂ 
                           
                             ( 
                             
                               
                                 I 
                                  
                                 
                                   ( 
                                   x 
                                   ) 
                                 
                               
                               - 
                               j 
                             
                             ) 
                           
                         
                          
                         
                             
                         
                          
                         
                           
                              
                             x 
                           
                           . 
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   10 
                   ) 
                 
               
             
           
         
       
     
         [0046]    With the method according to an exemplary embodiment of the present invention, local shape variations are found by deforming the multi-labeled image I T     N    through the transformation φ. This formulation allows an arbitrary number of regions to be segmented by optimizing only one function φ, in contrast to contour evolution methods, such as level set, where N functions are required to model contours (e.g., a level set function modeling each contour of a region Ω). The increasing number of contours in level set methods quickly becomes a complex memory problem. This problem is bypasses by encrypting the information of the different regions in a single multi-label image I T     N   . In addition, the method, in accordance with an exemplary embodiment of the present invention, provides a consistent structural relationship between the different regions where one transformation φ is optimized. 
         [0047]    Algorithm 1 (show below) describes how to compute the gradient of the similarity measure ∇S(I,I T     N   ∘φ). For each region, create a temporary binary image I i  of the region Ω i  and compute the corresponding probability density function p i . The image I i  is used when computing the gradient descent of this particular region ∇(I i ∘φ)(log p i (I(x))). The image I i  is chosen to be binary to avoid bias between different regions. The global gradient of the similarity measure of different regions is then updated.
   Algorithm 1 Similarity Measure for segmentation   Require: I,I T     N   =first approximation of N regions, φ.   Ensure: The gradient of the similarity measure ∇S(I,I T     N   ∘φ).   1: for Each region i in Ω do   2: Create a temporary image I i  corresponding to the region Ω i .   3: Compute p i  for the region Ω i  (equation (10)).   4: Compute ∇S(I,I i ∘φ)=∇(I i ∘φ)(log p i (l(x))).   5: Update ∇S(I,I T     N   ∘φ)+=∇S(I,I i ∘φ).   6: end for   
 
         [0057]    A description of experimental results of the multi-label segmentation method according to an exemplary embodiment of the present invention will now be provided. 
         [0058]      FIGS. 1A-C  show results of the segmentation, in accordance with an exemplary embodiment of the present invention. Here, eight different regions: liver, gallbladder, right kidney, left kidney, aorta, vena, cava, spleen and the background, were segmented. Image (a) in  FIGS. 1A-C  represents a rough initialization of I T     N    (hereinafter also referred to as T T     N   ) and image (b) in  FIGS. 1A-C  is a result of the multi-segmentation method according to an exemplary embodiment of the present invention, applied to its corresponding image (a). 
         [0059]    In image (b) of  FIG. 1A , six of the segmented regions are marked with an “X”. In image (b) of  FIG. 1B , four of the segmented regions are marked with an “X”. In image (b) of  FIG. 1C , three of the segmented regions are marked with an “X”. The marked regions in image (b) of  FIGS. 1A-C  clearly illustrate that the multi-label segmentation correctly delineates the different organs in the abdomen without leaking or overestimation. 
         [0060]    The liver segmentation result was compared to a ground-truth using five metrics: volumetric overlap, relative absolute difference, average symmetric absolute surface distance, symmetric RMS surface distance and maximum symmetric absolute surface distance. These metrics were evaluated using by assigning a score as described, for example, in [van Ginneken, B., Heimann, T., Styner, M.: 3d segmentation in the clinic: A grand challenge. In: 3D Segmentation in the Clinic: A Grand Challenge, MICCAI 2007 (2007), pp. 7-15]. Table 1 (shown below) presents the segmentation results. 
         [0000]    
       
         
               
               
               
               
               
               
               
             
               
               
               
               
               
               
               
             
           
               
                 TABLE 1 
               
               
                   
               
             
             
               
                 Metric 
                 V [%] 
                 Score 
                 dv [%] 
                 Score 
                 d moy  [mm] 
                 Score 
               
               
                   
               
               
                 Liver 
                 11.34 
                 57 
                 1.95 
                 90 
                 1.5 
                 60 
               
               
                   
               
             
          
           
               
                 Metric 
                 drms [%] 
                 Score 
                 d max  [%] 
                 Score 
                 Score 
                 total 
               
               
                   
               
               
                 Liver 
                 3.4 
                 50 
                 27.3 
                 65 
                   
                 64 
               
               
                   
               
             
          
         
       
     
         [0061]      FIGS. 2A  and B are flowcharts that illustrate a method for multi-label segmentation according to an exemplary embodiment of the present invention. 
         [0062]    In  FIG. 2A , an image I, an image I S  and pre-segmented labels I T     N    are input ( 205 ). In this example, the image I is a CT image of a patient&#39;s abdomen. It is to be understood, however, that this image could be of virtually any part of the patient&#39;s anatomy. In addition, this image could be have been acquired by a variety of imaging modalities, one such exemplary modality being magnetic resonance (MR). In this example, the image I S  is a baseline image that corresponds to a patient&#39;s abdomen. It is to be understood that image I S  is not the same image as image I. Further, image I S  has corresponding pre-segmented labels I T     N   . The pre-segmented labels I T     N    are a good segmentation of certain organs in the abdomen of the image I S . The pre-segmented labels I T     N    are manually marked by a doctor, for example. 
         [0063]    After the images I and I S  are input, they are aligned ( 210 ). This is done by using the fluid-based technique described by equations 5 and 6 with an LCC similarity measure, for example. The result of this alignment is a mapping/transformation φ*. This mapping/transformation φ* is applied to I T     N    to get T T     N    ( 215 ). For example, the warping is applied by using tri-linear interpolation, e.g., I T     N   ∘φ*. Hereinafter, I T     N   ∘φ* may be referred to just as T T     N   . In other words, T T     N    is a rough initialization of the pre-segmented labels I T     N    for the image I. As already mentioned, an example of this rough initialization is shown in image (a) of  FIGS. 1A-C . 
         [0064]    Now the roughly-initialized (e.g., deformed) pre-segmented labels image T T     N    is aligned to the image I by maximizing the likelihood of intensity distributions ( 220 ). In other words, the pre-segmented labels image T T     N    is updated with a new mapping/transformation φ until a desired refined segmentation of the organs is achieved. This process will now be described. 
         [0065]    Using the image I and the roughly-initialized pre-segmented labels image T T     N   , ν k , which is a gradient of the similarity measure ∇S(I,I T     N   ∘φ) (e.g., eq. (9)), is computed ( 225 ). This step will be described in more detail hereinafter with reference to  FIG. 2B . The gradient ν k  is regularized ( 230 ) with Gaussian smoothing. A new mapping/transformation φ is computed by applying the regularized gradient to eq. (5) ( 235 ). This can be seen as an instance of Christensen et al.&#39;s fluid registration, discussed previously. The new mapping/transformation φ is used to update the roughly-initialized pre-segmented labels image ( 240 ), e.g., by computing T T     N   ∘φ. The sequence of steps (outlined in  220 ) is repeated until the cost function of the similarity measure stops decreasing, for example. As already mentioned, an example of the results of aligning the pre-segmented labels image T T     N    to the image I is shown in image (b) of  FIGS. 1A-C . 
         [0066]    The left-hand side of  FIG. 2B  illustrates the process of computing ν k  in step  225 . This process is done for every label i. An example of several labels that will undergo this process is shown by 1, 2, 3, 4 and 5 (including the background identified as a separate region) identified as T T     N    on the right-hand side of  FIG. 2B . Using the image I and the deformed pre-segmented labels image T T     N    from box  215  (the example of which is shown on the right-hand side of this figure), a temporary image I i ∘φ for the region Ω i  is created ( 225   a ). The temporary image being I 1  for label  1  (i.e., region Ω 1 ). Using equation (10), the intensity distribution function for the region Ω i  is computed ( 225   b ). The gradient of the similarity measure of the temporary image ∇S(I,I i ∘φ)=∇(I i ∘φ)log p i (I(x))) is computed ( 225   c ). The final gradient of the similarity measure, i.e., the final gradient image ∇S(I,I T     N   ∘φ)+=∇S(I,I i ∘φ), is updated by concatenating the final gradient image with the gradients of the current label. This process is then repeated for I 2  for label  2  (i.e., region Ω 2 , I 3  for label  3  (i.e., region Ω 3 ), I 4  for label  4  (i.e., region Ω 4 ) and I 5  for label  5  (i.e., region Ω 5 ). A example of the different regions and temporary images for each label is shown by the shaded labels  1 ,  2 ,  3 ,  4  and  5  in images I 1 ,I 2 , I 3 , I 4  and I 5 , of  FIG. 2B , respectively. 
         [0067]    A system in which exemplary embodiments of the present invention may be implemented will now be described with reference to  FIG. 3 . As shown in  FIG. 3 , the system includes a scanner  305 , a computer  315  and a display  310  connected over a wired or wireless network  320 . The scanner  305  may be an MR or CT scanner, for example. The computer  315  includes, inter alia, a central processing unit (CPU)  325 , a memory  330  and a multi-label segmentation module  335  that includes program code for executing methods in accordance with exemplary embodiments of the present invention. The display  310  is a computer screen, for example. 
         [0068]    It is understood that the present invention may be implemented in various forms of hardware, software, firmware, special purpose processors, or a combination thereof. In one embodiment, the present invention may be implemented in software as an application program tangibly embodied on a program storage device (e.g., magnetic floppy disk, RAM, CD ROM, DVD, ROM. and flash memory). The application program may be uploaded to, and executed by, a machine comprising any suitable architecture. 
         [0069]    It is also understood that because some of the constituent system components and method steps depicted in the accompanying figures may be implemented in software, the actual connections between the system components (or the process steps) may differ depending on the manner in which the present invention is programmed. Given the teachings of the present invention provided herein, one of ordinary skill in the art will be able to contemplate these and similar implementations or configurations of the present invention. 
         [0070]    It is further understood that the above description is only representative of illustrative embodiments. For convenience of the reader, the above description has focused on a representative sample of possible embodiments, a sample that is illustrative of the principles of the invention. The description has not attempted to exhaustively enumerate all possible variations. That alternative embodiments may not have been presented for a specific portion of the invention, or that further undescribed alternatives may be available for a portion, is not to be considered a disclaimer of those alternate embodiments. Other applications and embodiments can be implemented without departing from the spirit and scope of the present invention. 
         [0071]    It is therefore intended, that the invention not be limited to the specifically described embodiments, because numerous permutations and combinations of the above and implementations involving non-inventive substitutions for the above can be created, but the invention is to be defined in accordance with the claims that follow. It can be appreciated that many of those undescribed embodiments are within the literal scope of the following claims, and that others are equivalent.