Computerized Medical Imaging and Graphics 38 (2014) 725–734

Contents lists available at ScienceDirect

Computerized Medical Imaging and Graphics journal homepage: www.elsevier.com/locate/compmedimag

Multi-region labeling and segmentation using a graph topology prior and atlas information in brain images Saif Dawood Salman Al-Shaikhli ∗ , Michael Ying Yang ∗ , Bodo Rosenhahn Institute for Information Processing, Leibniz University Hanover, Appelstr. 9A, 30167 Hannover, Germany

a r t i c l e

i n f o

Article history: Received 30 September 2013 Received in revised form 12 April 2014 Accepted 13 June 2014 Keywords: Segmentation Labeling Multi-region Atlas information Topological graph Multi-level set Medical image

a b s t r a c t Medical image segmentation and anatomical structure labeling according to the types of the tissues are important for accurate diagnosis and therapy. In this paper, we propose a novel approach for multi-region labeling and segmentation, which is based on a topological graph prior and the topological information of an atlas, using a modified multi-level set energy minimization method in brain images. We consider a topological graph prior and atlas information to evolve the contour based on a topological relationship presented via a graph relation. This novel method is capable of segmenting adjacent objects with very close gray level in low resolution brain image that would be difficult to segment correctly using standard methods. The topological information of an atlas are transformed to the topological graph of a low resolution (noisy) brain image to obtain region labeling. We explain our algorithm and show the topological graph prior and label transformation techniques to explain how it gives precise multi-region segmentation and labeling. The proposed algorithm is capable of segmenting and labeling different regions in noisy or low resolution MRI brain images of different modalities. We compare our approaches with other state-of-the-art approaches for multi-region labeling and segmentation. © 2014 Elsevier Ltd. All rights reserved.

1. Introduction Multi-region image segmentation is a major task in medical imaging and it is important in diagnosis and therapy [1]. Due to poor resolution and weak contrast, image segmentation is difficult in the presence of noise and artifacts [2]. Many existing methods for segmentation are based on image intensity information, shape properties or shape priors [2–5]. Many researches addressed that the medical imaging systems like MRI, although it relatively provides high-resolution anatomical details but the identification of tissue information is limited by several factors like noise and image non-uniformity due to magnetic field inhomogeneities [1]. This gives difficulties of the brain tissue labeling [1,6]. Manual labeling of brain structures is achieved using a lot of information including image intensities, anatomical landmarks, position relative to neighboring brain structures and global position within the brain [7], which need long processing time. Therefore, the automatic labeling is necessary and to maintain it for brain tissue we consider a topological graph information and tissue labeling for the segmentation that give a precise knowledge about position, size and type of the brain tissue [8]. The topological

∗ Corresponding authors. Tel.: +49 511 762 5319. E-mail address: [email protected] (S.D.S. Al-Shaikhli). http://dx.doi.org/10.1016/j.compmedimag.2014.06.008 0895-6111/© 2014 Elsevier Ltd. All rights reserved.

graph prior gives an information about the organs in the brain images and the atlas gives useful information about the label of the brain regions. The labels transfer from the atlas to the target image after warping the atlas with the target image. Okada et al. [9] proposed multi-organ segmentation of the upper abdomen by finding the interrelations between the organs based on canonical correlation analysis. Suzuki et al. [10] proposed an atlas based multi-organ segmentation and detection of missing organ in abdominal CT images. Shimizu et al. [11] proposed simultaneous extraction of multiple organs from abdominal CT using abdominal cavity standardization process with feature database and atlas guided segmentation incorporating parameter estimation for organ segmentation. Linguraru et al. [12] proposed multi-region segmentation using graph cut method for four abdominal organ segmentation. Kohlberger et al. [13] proposed multi-organ segmentation from CT medical images using learning-based segmentation and shape representation. Bazin and Pham [14] proposed multiregion segmentation algorithm of brain image using topological and statistical atlases of brain as prior to the segmentation framework. Nocera and Gee [7] proposed tissue classification of cerebral magnetic resonance images using Bayesian estimation method. Fischl et al. [15] proposed an automatic labeling of Neuroanatomical structures in the human brain by estimating the probability information from manual labeled training data. Soni [16] proposed brain tissue classification of only three types of tissue (gray matter, white

726

S.D.S. Al-Shaikhli et al. / Computerized Medical Imaging and Graphics 38 (2014) 725–734

Fig. 1. Block diagram explain the proposed algorithm for multi-region labeling and segmentation.

matter and CSF) using conditional random field for magnetic resonance images. Cocosco et al. [17] proposed a full automatic brain tissue classification method for three types of brain tissue in magnetic resonance images by measuring a tissue probability map. Liu et al. [18] proposed a method for image segmentation using multi-context label tree structure. Sabuncu et al. [19] proposed a nonparametric image segmentation and label fusion approach using image registration. Aljabar et al. [20] proposed a framework for brain tissues multi-atlas segmentation but the accuracy of this approach degrades in the presence of strong noise and by using single atlas information. Mansouri et al. [21] proposed multi-region competition algorithm for intensity-based image segmentation. Vazquez et al. [22] proposed image segmentation algorithm from the viewpoint of image data regularized clustering. Aforementioned work of multi-region segmentation and labeling focused on either image region intensity with atlas topology [14] or region priors [21] or multi-atlas information for segmentation [20], however the shape of the region of brain tissues may vary from person to person and the intensity differs according to the image modality. Also a multi-atlas methods need more training data and more computational time. Furthermore, the performance of these works degrade in the presence of high level of noise. In contrast to the aforementioned works, our contribution of this paper is multi-region segmentation and labeling using a multilevel set formulation which includes a topological graph prior and atlas information in brain images. We propose a new topological representation of brain regions based on object-oriented graph to define an image and use this topological representation to constrain the level set functional energy. Also the proposed object-oriented graph representation is used to compute the similarity between

target and atlas template images to map the topological information from an atlas to a target image. The topological graph prior is embedded in the multi-level set energy equation and acts as an additional prior term to identify both the overlapped regions and weak boundaries between adjacent regions in the image, as shown in Fig. 1. The topological graph prior allow us to handle the huge variability of different modalities of the brain images. Consequently, our algorithm is less sensitive to noise and gives accurate segmentation of ambiguous regions of the brain. For brain segmentation and labeling, we propose seven labels of brain tissue as shown in Fig. 3 and Table 3. The outcome of our algorithm is conjoint of multi-class image segmentation and labeling. In all of our experiments, we concentrate on brain segmentation, however, it is worth noting that the method is general and can be applied to other scenarios, for example abdominal organ segmentation by computing the topological relationship of the abdominal organs using abdominal atlas. The organization of this paper is as follows: Section 2 explains the proposed approach. The discussion of the experimental results is presented in Section 3. Finally, the conclusion of the paper is presented in Section 4.

2. Method In this section we will explain our proposed method for multiregion segmentation and labeling based on a multi-level set formulation with an atlas information and a topological graph prior. The proposed method is based on three steps; Topological graph prior step, label transformation step, and level set energy minimization step.

S.D.S. Al-Shaikhli et al. / Computerized Medical Imaging and Graphics 38 (2014) 725–734

727

Fig. 2. Example of image representation as topological graph.

2.1. Graph prior

a22 as primary conditions of the topological relationship. This is summarized as follows:

Human body organs have specific topological correlations between them and according to these correlations, the exact location and boundary of these organs can be determined. If we consider the image B as sets of clusters (segments) B = Oi , Oi+1 , . . ., ON depending on the dissimilarity between them and Oi is the membership function of each cluster. These clusters are connected with each other by a specific topological relationship then the topological graph of these clusters can give information like the area, the location and the topological relationship of each cluster in the image. The topological graph is constructed from test image to provide the prior knowledge to the segmentation process. The clusters in the topological graph of the image B are determined using Otsu’s method [23] and these clusters are labeled according to their topological relationship. We consider four types of topological relations (disjoint, contact, inside, overlap). Egenhofer and Herring [24] have defined and computed O◦ as an interior of the cluster, Oc as a complement (exterior) of the cluster and ∂O as a boundary of the cluster. The topological region relationship (TRL) between the clusters is calculated in terms of probability of intersections of these clusters in a 9-intersection model in 3 × 3 matrix, as follows:

⎧ a11 = 0, a12 = 0, a21 = 0, a22 = 0 if RLdis (Oi , Oi+1 ) > 0 ⎪ ⎪ ⎪ ⎪ ⎨ a11 = 0, a12 = 0, a21 = 0, a22 = 1 if RLcon (Oi , Oi+1 ) > 0 ⎪ a11 = 0, a12 = 0, a21 = 1, a22 = 0 if RLin (Oi , Oi+1 ) > 0 ⎪ ⎪ ⎪ ⎩

a11 = 1, a12 = 1, a21 = 1, a22 = 1 if RLov (Oi , Oi+1 ) > 0

where RLdis , RLcon , RLin , and RLov are disjoint, contact, inside, and overlap region relationship respectively, as follows: RLdis (Oi , Oi+1 ) = 1 − max{|Oi (b) + Oi+1 (b) − 1|} RLcon (Oi , Oi+1 ) = min{(1 − max(|O◦ (b) + O◦ (b) − 1|)), b

max(min(∂O (b), ∂O i

b



◦

a11 Oi

a31 Oi



◦

Oi+1

Oi+1





◦

a12 Oi

a32 Oi



∂Oi+1

∂Oi+1





◦

a13 Oi

a33 Oi



c Oi+1

Oi+1

(1)

Each element in Eq. (1) represents specific topological relationship. For example, to achieve the disjoint relation, we need to a11 = 0, a12 = 0, a21 = 0 and a22 = 0, which means that all pixels in cluster Oi are not in Oi+1 . To achieve the inside relation between Oi and Oi+1 , we need to a21 = 1, which means that the pixels in ◦ ∂Oi are in Oi+1 . To achieve the contact relation, we need to a11 = 0, a22 = 1, a12 = 0, and a21 = 0 which mean that the pixels in ∂Oi are in ∂Oi+1 . The overlapped region is achieved by a11 = 1, a12 = 1, a21 = 1, and a22 = 1. Therefore, we consider the elements a11 , a12 , a21 , and

i+1

i

i+1

(b)))}

(4)

RLin (Oi , Oi+1 ) = min(1, min(1 + O◦ (b) − Oi (b))) b

(5)

i+1

RLov (Oi , Oi+1 ) = min{max(min(O◦ (b), O◦ (b))), b

b

 ⎞      c ⎟ ⎜   ◦  a22 ∂Oi ∂Oi+1 a23 ∂Oi Oi+1 = ⎝ a21 ∂Oi Oi+1 ⎠  c ◦   c   c c  ⎛

(3)

b

max(min( TRL(Oi , Oi+1 )

(2)

◦ O i

(b), ∂O

i+1

max(min(∂O (b), ∂O b

i

i+1

i+1

i

(b))), max(min(O◦ (b), ∂O (b))), b

(b)))}

i+1

i

(6)

Table 1 summarizes how each element of the matrix in Eq. (2) determines the relationship between the clusters by checking the primary conditions and the secondary conditions. The primary conditions are the main conditions to determine the topological relationship between the regions. The secondary conditions are proposed to be ‘ones’ in the matrix. The secondary conditions in Eq. (1) represent the intersection of the region of interest (ROI) with the complement of the other region, which means the intersection of the ROI with itself. For example, Oi and Oi+1 are two regions in the image. If we choose the element a13 in Eq. (1) as an example of the secondary conditions, which represents the intersection of Oi with the complement c . O will intersect with itself because O is part of Oc . Thereof Oi+1 i i i+1 fore always the intersection is achieved in the secondary conditions and we set them to ‘ones’ in the matrix.

728

S.D.S. Al-Shaikhli et al. / Computerized Medical Imaging and Graphics 38 (2014) 725–734

Table 1 Topological properties of different regions of the image according to Eqs. (1) and (2). TRL(Oi , Oi+1 )

Primary conditions

Secondary conditions

Contact regions Inside regions Disjoint regions Overlap regions

a22 = 1, a11 = 0, a12 = 0, a21 = 0 a21 = 1, a11 = 0, a12 = 0, a22 = 0 a11 = 0, a12 = 0, a21 = 0, a22 = 0 a11 = 1, a12 = 1, a21 = 1, a22 = 1

a13 = 1, a23 = 1, a31 = 1, a32 = 1, a33 = 1 a13 = 1, a23 = 1, a31 = 1, a32 = 1, a33 = 1 a13 = 1, a23 = 1, a31 = 1, a32 = 1, a33 = 1 a13 = 1, a23 = 1, a31 = 1, a32 = 1, a33 = 1

The topological similarity Ts between each cluster in the topological graph and the corresponding region in the image during evolution is for updating the labels of each region in the image at each t during evolution process, i.e. it is applied iteratively during curve evolution to update the label of each region in an image. Ts is determined by subtraction TRLtotal (Oi , ROi ) from TRLtotal (Ri , RRi ) during evolution process:

 Ts (Oi , Ri ) =

Fig. 3. Topological graph of the atlas template with its labels.

Fig. 2 explains the representation of the anatomical structures in the image as topological graph. The red rectangle in Fig. 2 shows that the labels that are considered in our calculation depend on the atlas labels illustrated in Fig. 3 (as we will explain in the next subsection), in another words, the labels may consist of any combination depending on the slice level of the brain image. Table 2 shows TRL of each region in the image of Fig. 2. The connected components represent the total relationship of each region. The number of cavities in each region indicates how many regions are inside it or held by it. For example, in case of region white matter (WM) in Fig. 2, it has two relationships (two connected components), (1) TRL1 (WM, GM) and (2) TRL2 (WM, ventricles) as in Eq. (7). TRL1 (WM, GM) shows that WM is inside GM, while TRL2 (WM, ventricles) shows that WM covers ventricles according to Eq. (1):

⎛

0

TRL1 (WM, GM) = ⎝ 1 1

0 0

1

⎞

0

TRL2 (WM, ventricles) = ⎝ 0 1

1 1 0

⎞

1⎠

(7)

1 1

TRLtotal (WM, RWM ) = TRL1 (WM, GM) + TRL2 (WM, ventricles)

⎛

=

0

1

2

2 2

2

R

i

)

(9)

otherwise

where (Ri ) are the regions inside their contours during evolution process. ROi and RRi are the regions that have topological relationship with Oi and Ri , respectively. Also the area and the centroid of each contour are calculated at each t during evolution process and compared with the area and the centroid of corresponding cluster in the topological graph of the image:

 Ai =

dA,

1 A

Cxi =

A

 xe dA,

Cyi =

A

1 A

 ye dA

(10)

A

where Ai are the areas, Cxi and Cyi are the coordinates of centroid, xe and ye are coordinates of the centroid of the differential element of area dA. These prior information are added to the functional energy as topological graph prior term:

 N−1

Eg [(i )i=1 ] = ˛(

(|AOi − ARi |)dx + (|COi − CRi |) + Ts ) R

i



(11)



Topological graph prior term

1 1

⎛

1



1⎠,

for TRLtotal (Oi , ROi ) = TRLtotal (Ri , R

0

⎞

⎝1 0 2⎠

(8)

where TRLtotal is the total topological relationship, RWM are the regions that have topological relationship with (WM). In Eq. (8), a12 = 1 indicates that the region WM covers ventricles while a21 = 1 indicates that WM is inside GM. Table 2 Topological properties of different regions of the image in Fig. 2. Region label

Region name

# of connected components/(region name)

Internal cavity/handles

1 2 3 4 5

Skull Sulcal CSF GM WM Ventricles

2/(bg), (sulcal CSF) 2/(skull), (WM) 2/(sulcal CSF), (WM) 2/(GM), (Ventricles) 2/(sulcal CSF), (WM)

1/(Sulcal CSF) 1/(GM) 1/(WM) 1/(Ventricles) 0

Eg is the energy of the topological graph. ˛ is constant (˛ = 1 or 0) to run the algorithm with or without topological graph prior. AOi , COi are the area and centroid of the clusters in topological graph and ARi , CRi are the area and the centroid of the regions in the image B during the evolution process respectively. As mentioned above, the accuracy of the segmentation depends on the accuracy of the extraction of the topological graph information which may be affected in the presence of high level of noise. Therefore we propose to use an atlas information as an additional prior information to solve this problem, as we will explain in the next sections. 2.2. Construction of atlas template and its topological properties The segmentation accuracy depends on the accuracy of the topological graph prior extraction. In the presence of strong noise, the extraction of the topological graph will be affected. Therefore, we propose to use an atlas information to increase the accuracy of the topological graph in presence of strong noise. The atlas template T is constructed using the average of 30 brain dataset from brainweb [25]. The average atlas template is calculated using average brain model [26]. The topological properties of the atlas are summarized in Fig. 3 and Table 3. Fig. 3 shows the topological graph of the atlas template with the label of each region and Table 3 explains the topological properties of the atlas template.

S.D.S. Al-Shaikhli et al. / Computerized Medical Imaging and Graphics 38 (2014) 725–734 Table 3 Topological properties of the atlas template.

2.5. Multi-level set formulation and curve evolution

Tissue label

Tissue type

# of connected components

Internal cavity/handles

1 2 3

Sulcal CSF Cortical gray matter (GM) White matter (WM)

3 2 3

4 5 6 7

Subcortical gray matter Ventricles Cerebellum Brain stem

2 3 2 1

1/(GM) 1/(WM) 2/(Ventricles), (Subcortical GM) 0 0 0 0

This section describes the proposed (modified) multi-level set N method with labeled topological graph prior. Let B = ( i=1 Ri ) is the input image with N regions, i ∈ [1, 2, . . ., N]. We assume that for each Ri there is its complement Ric : Ri (t) = {u ∈ R|i (u, t) > 0}, c



(12)

where O are the objects in T. Assuming that the topological properties of the atlas template are summarized in Fig. 3 and Table 3, while the topological properties of B are calculated as explained in Section 2.1. After constructing the topological graph for both the atlas template (atlas training data) and the input image, the similarity between T and B with associated graphs G1 (OB , EB ) and G2 (OT , ET ) is measured. OB are the objects in B, OT are the objects in T. EB and ET are the edges between the objects in both B and T, respectively. The object similarity and the topological similarity are considered to determine the overall similarity between B and T as follows:



On

(wOn wOi Oe (On , Oi ))

w ∈T On



(13)

w O ∈B Oi i

Oe (On , Oi ) = 1 − max{|On (u(x, y)) − Oi (u(x, y))|}

where Os is the object similarity, n = 1, 2, . . ., 7, wOn and wOi indicate the importance given to On and Oi while computing the similarity. Oe is the object equality. u(x, y) ∈ On and Oi . The topological similarity between B and T is determined by measuring the similarity of TRLtotal (On , ROn ) and TRLtotal (Oi , ROn ). To find the best matching between the topological graph of B and the atlas template T, the similarity function SIM is measured: SIM(B, T ) = Os (On , Oi ) + Ts (On , Oi )

(15)

2.4. Registration of the atlas information and label transformation



(|AOl − ARi |)dx + (|COl − CRi |) + Ts ) i

i

i



(18)

Ri (t)



ωi (b)db +

i (b)db +  Rc

R

i

i

i

ds + Egl

(19)

where the first two terms are the data terms of the region. The third term is the regularization term and the fourth term is our proposed prior term.  is positive real constant to weigh the relative contribution of the energy equation. ωi are the data in Ri and i are the data in Ric . The data term is modified to be constrained by Ts . is con2 strained by Ts : ωi (b) = [− ln Pi (B(b)) = (B(b) − i ) ] and i (b) = 2

[− ln Pj (B(b)) = (B(b) − j ) ]Ts where i are the mean over Ri and j are the mean over Rc . According to Eqs. (9) and (19), Ts i

represents the topological function of the label state of the set Ric . To minimize Eq. (11) by curve evolution we compute: di ∂E =− dt ∂i

(20)

∂E are the derivative of functional energy with respect to  and i ∂i

they are computed as for the standard region computation functional in [28]. Following [28], we get the evolution equation of the curves i :

⎛

⎞

⎜ ∂i = − ⎝ωi (b) − ∂t

i (b)



 



⎟

+ ˛[ |AOl − ARi | + |COl − CRi | + Ts ] + ki ⎠ ni



i



i



(21)

Labeled topological graph prior

where ki are the curvature of zero level set of i , ni are the external unit normal of the curve, i ∈ [1, . . ., N], j ∈ [1, . . ., N] and i = / j. During curve evolution, the curves are constrained by the labeled topological graph prior term and the curvature term. For N-region segmentation, let Ri be regions in the image B (Ri ∈ B) and let i (0) be an initial curve and i (t) is a curve in an iteration t: i



R

i=1

(17)

1. AR and CR are updated for each time step during evolution

The human organs are non-rigid organs and normally there is a relative shape difference from person to person. Sometimes due to high level of noise or low image resolution the topological graph gives limited information about the regions in the image. Therefore a multi-modality non-rigid demon algorithm is proposed to use for image registration and label transformation from the atlas to the target image [27]. After label transformation from the atlas to the topological graph of the target image, Eq. (11) should be rewritten in a form of label topological graph prior: N



N Etotal [(i )i=1 ] =

(14)

u

Egl [(i )i=1 ] = ˛(

c

c 

 N

The total Euler–Lagrange energy functional can be written as follows:

T = O1 , O2 , . . ., O7

On ∈T, Oi ∈B

c

R1 (t), R1 (t) ∩ R2 (t), R1 (t) ∩ R2 (t) ∩ R3 (t), . . .,

We consider the graph similarity (object similarity and topological relationship similarity) as a similarity measurement between the test image and the atlas template. The topological graph in T are constructed by measuring the objects or regions in each image in atlas template using Otsu’s method [23] and then determine the topological relationship between these regions:

Os (On , Oi ) =

i = 1, . . ., N



2.3. Selection of appropriate atlas template



729

(16)



Labeled topological graph prior term

where l is the label of the region as explained in Fig. 3 and Table 3.

i

process and compared with the area and the centroid of the labeled topological graph prior. The errors between AOi and AR and between COi and CR should be minimized. i

i

2. The topological similarity Ts defines the label state of each region in the image B at time t + 1 with respect to the label of the same region at time t. Using Eq. (9), if Ts = 0 at i (t) and i (t + 1) then the label of the region is still without changing, i.e. Ri ∈ i . If Ts = 0 at i (t) and Ts = 1 at i (t + 1) then the label of the region will change, i.e. Rj ∈ j , i = / j, i ∈ [1, . . ., N] and j ∈ [1, . . ., N]. 3. If b is a point of contact between two curves (i , j ), then the curve will be constrained by the curvature term as follows: If the curvatures are positive (ki (b) ≥ 0, kj (b) ≥ 0) this indicates that these curves are retract and not intersect. If (ki (b) ≤ 0, kj (b) ≥ 0) this indicates that these two curves will be in the same direction but because |ki (b) ≤ kj (b)|, the curve j retracts faster than i and the curves will not intersect. The algorithm of the topological graph prior is described in Algorithm 1. The graph constraint

730

S.D.S. Al-Shaikhli et al. / Computerized Medical Imaging and Graphics 38 (2014) 725–734

makes the partitioning more precise during evolution process by adding addition constrain information (Ts , area, centroid). Algorithm 1. graph prior

Algorithm for computing the labeled topological

Given an image B, consists of N regions. 1. Label transformation from the atlas to the topological graph after non-rigid registration [27] 2. Compute the error of A and C between Ri and Oi during evolution process: • During curve evolution at each t. While (|AOl − AR | > )&(|COl − CR | > ) do i

i

i

˛[(|AOi − AR |) + (|COi − CR |)]ni l

i

l

i

Table 4 Segmentation accuracy for each database without the effect of noise. Algorithm

DSC [29]

DSC [30]

DSC [25]

Topological graph prior with atlas information (proposed) Topological graph prior without atlas information [8] With atlas information without topological graph [14] Atlas based segmentation [20] Without graph prior [21,22] ChanVese [3] Level set fuzzy based [5]

94%

91.9%

95.5%

93.56%

90.57%

94.88%

88.5%

89.5%

92%

90% 80.64% 61.82% 40.87%

87% 79.89% 61.6% 40.1%

93% 82.49% 62.78% 42.63%

i

end while 3. Compute the similarity of the topological relationship Ts and label update: • During curve evolution (at each t) for each region R ∈ B. • At time = t: ∀Ri ∈ B; if Ri > 0 then R ∈ i ; Find TRLtotal (Ri ); [li ] ={b(x, y) ∈ Ri |b ∼ l} ; end if • At time = t+1: if (TRLtotal (Oi , ROi ) − TRLtotal (Ri , RRi ) = 0;) then Ts = 0 ;  comment: Ri > 0 R ∈ i lt+1 = lt ; [li ] ={b(x, y) ∈ Ri |b ∼ l} ; / 0;) if (TRLtotal (Oi , ROi ) − TRLtotal (Ri , RRi ) = Ts = 1 ; / i ;  comment: Ri < 0 R∈ R ∈ j ;  comment: Rj > 0 / lt ; lt+1 = [lj ] ={b(x, y) ∈ Rj |b ∼ l} ; i=1:N;j=1:N;i = / j;

end if end if

3. Experiments To highlight the advantages of the proposed algorithm, the algorithm is tested in the presence of strong noise and with low resolution images for more than 300 images using the MedPix [29], Wesky E Snyder [30], the brain web for simulated brain database [25] and other medical images from the internet. Our algorithm is compared to state-of-the-art methods. The medical images in all databases are 2D MRI images and CT images. The sizes of the images are 181 × 217 and 512 × 512 for brain sections of MRI images and CT images. The experiments run in MATLAB using a 2.0 GHz Intel core I3 CPU. First, we will show the qualitative results of topological graph prior without atlas information. Then we will show the qualitative results of topological graph prior with atlas information. Finally we will show the quantitative evaluation of our algorithm. 3.1. Qualitative results of topological graph prior without atlas information The results of the first part of our algorithm show the improvement in the segmentation by using the topological prior information [8]. This information help the contour, during the evolution process, to detect the ambiguous regions in the image. Fig. 4 shows multi-region segmentation for abdominal and brain MRI images with and without topological graph prior. The ground truth is obtained by manual segmentation. Fig. 4 also shows the improvements of our algorithm to capture the overlapped and close gray level regions according to its topological location in the image. The abdominal image in Fig. 4 shows the improvements of our algorithm mainly in the segmentation of aorta, liver and diaphragm. The brain images show the segmentation of the cerebellum, brainstem, white matter and gray matters. In Fig. 4, the segmented regions

The significant of the bold values is to highlight the highest accuracy value.

in the proposed algorithm are labeled by colors according to their topological relationship. The result of the algorithms proposed by [21,22] are labeled manually to visualize the differences. The accuracy of this part of our algorithm depends on the accuracy of the precise extraction of each cluster in the topological graph, i.e. Ts , A and C should be computed precisely for each cluster in the topological graph. 3.2. Qualitative results of topological graph prior with atlas information In this section, we will show the result of using the atlas information as an additional prior information with the topological graph prior. As mentioned previously, the topological graph prior may affected in the presence of high level of noise. This part of the proposed algorithm solves this limitation. We propose to use an atlas template to label the graph of the input image and eliminate the effect of noise. Fig. 5 shows the results of using the atlas information with the topological graph compared with the results explained in Section 3.1 for multi-region segmentation [8] as well as with the approaches proposed in [14,20–22]. In Fig. 5 we observe the improvement of the segmentation in the presence of noise. The atlas information provides an accurate extraction of the topological graph information which improves the multi-region segmentation and labeling in noisy and low resolution images with less computational time compared to the state-of-the-art methods as illustrated in Tables 4–6. Fig. 5 demonstrates the segmentation results in the presence of high level of noise; the black rectangle (a) shows the input noisy images and the ground truth segmentation, the red rectangle (b) shows the segmentation with topological graph prior and atlas information with accurate region segmentation and labeling. The green rectangle (c) shows the segmentation with topological graph prior. The blue rectangle (d) shows the result of atlas registration based segmentation [14], and the yellow rectangle (e) shows the segmentation without any prior information [21,22]. The improvement of the segmentation of different brain tissues like white and gray matter, ventricles, and cerebellum can be seen. In Fig. 5, the skull is not signed in the proposed labels and it is segmented and labeled randomly according topological relationship with the other regions in the image. 3.3. Quantitative evaluation To validate the accuracy of our algorithm, we compare our algorithm with other state-of-the-art methods [3,5,8,14,20–22] using dice similarity coefficients (DSC) [31]. DSC is measured by computing the similarity between the ground truth segmentation and our algorithm as well as the methods proposed in [3,5,8,14,20–22]. A large DSC indicates higher accuracy: DSC(Igt − It ) =

2O(Igt − Itest ) O(Igt ) + O(Itest )

(22)

S.D.S. Al-Shaikhli et al. / Computerized Medical Imaging and Graphics 38 (2014) 725–734

731

Fig. 4. Multi-region segmentation results. (a) Input images, (b) ground truth, (c, d) proposed algorithm with graph prior and (e, f) without graph prior.  = 0.2, iteration = 70. An example images from database [29,30,25].

Table 5 Overall segmentation accuracy without the effect of noise of all images in database [29,30,25] and the average computation time for each frame. Algorithm

DSC

MAD

HD

# Iteration

Time

Segmentation with topological graph prior and atlas information (proposed) Segmentation with topological graph prior and without atlas information [8] Segmentation with atlas information and without topological graph [14] Atlas based segmentation [20] Segmentation without graph prior [21,22] ChanVese [3] Fuzzy based level set segmentation [5]

93.8% 93% 90% 90% 81% 62% 41.2%

0.66 mm 0.68 mm 0.90 mm 0.93 mm 1.9 mm 3.2 mm 5.3 mm

2.9 mm 3.2 mm 3.9 mm 4.5 mm 5.5 mm 7.2 mm 9.7 mm

70 70 70 70 70 400 300

2.45 min 2.24 min 3.1 min 2.9 min 2.87 min 5.15 min 3.2 min

where O(Igt − Itest ) is the number of overlapping pixels, O(Igt ) + O(Itest ) is the summation of the number of pixel in each image. We also employ the symmetric mean absolute distance (MAD) and Hausdorff distance (HD) [32] between the resulting segmentation and the corresponding reference segmentation as additional metrics to evaluate the segmentation results. MAD is calculated by measuring the average distance from all points on the border of the automatically segmented brain tissue to the border of the reference segmentation. On the other hand, to assess the maximal local discrepancy between an automatic segmentation and reference segmentation, the symmetric Hausdorff distance between the

border of the automatically segmented brain tissue and that of the reference segmentation is calculated. The smaller the MAD or Hausdorff distance, the better aligned the points on the two border and thus the better the agreement with the reference segmentation. Our algorithm is robust with respect to the level of noise and the number of the segmented region. Fig. 6 shows the stability of our algorithm as the number of the segmented region increases comparing with the other methods [3,5,8,14,20–22]. The proposed algorithm is based on a single topological atlas with the graph prior and it outperforms the other multi-atlas based segmentation [14,20]. The methods in [14,20] need more training data and more computational time because if every atlas is registered with the

Table 6 Overall segmentation accuracy with the effect of noise (standard deviation 0.16) of all images in database [29,30,25] and the average computation time for each frame. Algorithm

DSC

MAD

HD

# iteration

Time

Segmentation with topological graph prior and atlas information (proposed) Segmentation with topological graph prior and without atlas information [8] Segmentation with atlas information and without topological graph [14] Atlas based segmentation [20] Segmentation without prior information [21,22] ChanVese [3] Fuzzy based level set segmentation [5]

83% 70% 73% 77% 49% 20% 17.7%

1.0 mm 1.6 mm 2.5 mm 2.3 mm 5.1 mm 7.3 mm 7.5 mm

4.3 mm 5.0 mm 5.2 mm 5.1 mm 6.8 mm 9.9 mm 9.8 mm

70 70 70 70 70 400 300

2.45 min 2.24 min 3.1 min 2.9 min 2.87 min 5.15 min 3.2 min

732

S.D.S. Al-Shaikhli et al. / Computerized Medical Imaging and Graphics 38 (2014) 725–734

Fig. 5. Examples of a multi-region labeling and segmentation using database [25,29], (a) are the input images with noise and the ground truth, (b) are the segmentation with topological graph and atlas information, (c) are the segmentation with topological graph without atlas information [8], (d) are the segmentation with atlas information without topological graph [14], (e) are the segmentation without any prior [21,22],  = 0.2, iteration = 70, noise (SD = 0.16).

S.D.S. Al-Shaikhli et al. / Computerized Medical Imaging and Graphics 38 (2014) 725–734

733

Fig. 6. Effect of Gaussian noise on segmentation performance of database [25,29,30].

target image, the computational time of segmentation increases linearly with the size of the training data. In Fig. 6(a), we can see the performance of the proposed algorithm and the algorithms [14,20] is quite similar for 3-region segmentation but with increase of the level of noise the performance of the proposed algorithm is more robust than in [14,20]. Also as the number of the segmented region increases the performance of the algorithm [14,20] decreases comparing with our algorithm as explained in Fig. 6(b). Table 4 shows the accuracy (DSC) of our algorithm compared to the other methods for each database [25,29,30]. Tables 5 and 6 show the overall accuracy (DSC, MAD, and HD) of our algorithm and other algorithms with and without presence of noise using the images in the databases [25,29,30]. For all these investigated scenarios, our algorithm outperforms other methods. 4. Conclusion and future work We propose to use a topological graph prior with atlas information in a modified multi-level set formulation for multi-region segmentation and partitioning in brain images. As a high-level prior, it gives accurate region partitioning with respect to their topological location and relationship as well as with the atlas information. The accuracy of the proposed approach depends on the

accuracy of the extraction of the topological graph prior information which is achieved using the atlas information. The proposed method has higher accuracy while less computation time than state-of-the-art methods. Acknowledgement The work was partially funded by DAAD scholarship (A/10/96106) and MOHESR-Iraq (Baghdad University). The authors gratefully acknowledge these supports. References [1] Shattuck DW, Sandor-Leahy SR, Schaper KA, Rottenberg DA, Leahy RM. Magnetic resonance image tissue classification using a partial volume model. Neuroimage 2001;13(1):856–76. [2] Andrews S, McIntosh C, Hamarneh G. Convex multiregion probabilistic segmentation with shape prior in isometric log-ratio transformation space. In: ICCV. 2011. p. 2096–103. [3] Chan T, Vese L. Active contours without edges. IEEE Trans Image Process 2001;10(2):266–77. [4] Rathke F, Schmidt S, Schnoerr C. Order preserving and shape prior constrained intra-retinal layer segmentation in optical coherence tomography. In: MICCAI 2011. 2011. p. 370–7. [5] Li BN, Chui CK, Chang S, Ong SH. Integrating spatial fuzzy clustering with level set methods for automated medical image segmentation. Comput Biol Med 2011;41:1–10.

734

S.D.S. Al-Shaikhli et al. / Computerized Medical Imaging and Graphics 38 (2014) 725–734

[6] Kapur T, Grimson WEL, Wells WM, Kikinis R. Segmentation of brain tissue from magnetic resonance images. Med Image Anal 1996;1:109–27. [7] Nocera L, Gee JC. Robust partial volume tissue classification of cerebral MRI scans. In: Hanson KM, editor. SPIE Medical Imaging, vol. 3034. 1997. p. 312–22. [8] Al-Shaikhli SDS, Yang MY, Rosenhahn B. Medical image segmentation using multi-level set partitioning with topological graph prior. In: Pacific-Rim Symposium on Image and Video Technology PSIVT Workshops, LNCS, vol. 8334. 2013. p. 157–68. [9] Okada T, Linguraru MG, Yoshida Y, Hori M, Summers RM, Chen YW, et al. Abdominal multi-organ segmentation of CT images based on hierarchical spatial modeling of organ interrelations. In: MICCAI 2011. 2011. p. 173–80. [10] Suzuki M, Linguraru MG, Summers RM, Okada K. Analyses of missing organs in abdominal multi-organ segmentation. Abdom Imaging Comput Clin Appl 2012;7029:256–63. [11] Shimizu A, Ohnoa R, Ikegamia T, Kobatakea H, Nawanob S, Smutekc D. Segmentation of multiple organs in non-contrast 3D abdominal CT images. Int J Comput Assist Radiol Surg 2007;2:135–42. [12] Linguraru MG, Pura JA, Pamulapati V, Summers RM. Statistical 4D graphs for multi-organ abdominal segmentation from multiphase CT. In: MICCAI 2010. 2010. p. 89–96. [13] Kohlberger T, Sofka M, Zhang J, Birkbeck N, Wetzl J, Kaftan J, et al. Automatic multi-organ segmentation using learning-based segmentation and level set optimization. In: MICCAI 2011, vol. 6893, no. 3. 2011. p. 338–45. [14] Bazin PL, Pham DL. Homeomorphic brain image segmentation with topological and statistical atlases. MICCAI 2007. Med Image Anal 2008;12(5):616–25. [15] Fischl B, Salat DH, Busa E, Albert M, Dieterich M, Haselgrove C, et al. Whole brain segmentation: automated labeling of neuroanatomical structures in the human brain. Neuron 2002;33:341–55. [16] Soni A. Brain tissue classification of magnetic resonance images using conditional random fields. Department of Computer Sciences, University of Wisconsin-Madison; 2007, May. [17] Cocosco CA, Zijdenbos AP, Evans AC. A fully automatic and robust brain MRI tissue classification method. Med Image Anal 2003;7:513–27.

[18] Liu L, Zhang Y, Zhang C. Medical image segmentation based on contextual label tree for retrieval. J Comput Inform Syst 2011;7:1472–8. [19] Sabuncu MR, Yeo BTT, Leemput KV, Fischl B, Golland P. A generative model for image segmentation based on label fusion. IEEE Trans Med Imaging 2010;29:1714–29. [20] Aljabar P, Heckemann RA, Hammers A, Hajnal JV, Rueckert D. Multi-atlas based segmentation of brain images: atlas selection and its effect on accuracy. Neuroimage 2009;46:726–38. [21] Mansouri R, Mitiche A, Vazquez C. Multiregion competition: a level set extension of region competition to multiple region image partitioning. Comput Vis Image Underst 2006;101:137–50. [22] Vazquez C, Mitiche A, Ayed IB. Image segmentation as regularized clustering: a fully global curve evaluation method. In: ICIP 2004. 2004. p. 3467–70. [23] Otsu N. A threshold selection method from gray-level histograms. IEEE Trans Syst Man Cybern 1979;9(1):62–6. [24] Egenhofer M, Herring J. Categorizing binary topological relations between regions, lines and points in geographic databases. Technical report. Dept. of Surveying Eng., Univ. of Maine; 1991. [25] Cocosco CA, Kollokian V, Kwan RK-S, Evans AC. BrainWeb: online interface to a 3D MRI simulated brain database. Neuroimage 1997;5(4):S425. [26] Guimond A, Meunier J, Thirion JP. Average brain models: a convergence study. Comput Vis Image Underst 2000;77(2):192–210. [27] Kroon DJ, Slump CH. MRI modality transformation in demon registration. In: ISBI 2009, IEEE Proceeding. 2009. p. 963–6. [28] Zhu SC, yuille A. Region competition: unifying snakes, region growing, and Bayes/mdl for multiband image segmentation. PAMI 1996;18:884–900. [29] MedPix(r) Medical Image Database; 1999. [30] Wesley E Snyder: NC State University Image Analysis Laboratory Database; 2002. [31] Zou KH, Warfield SK, Bharatha A, Tempany CMC, Kaus MR, Haker SJ, et al. Statistical validation of image segmentation quality based on a spatial overlap index. Acad Radiol 2004;11:178–89. [32] Wang J, Li F, Li Q. Automated segmentation of lungs with severe interstitial lung disease in CT. Med Phys 2009;36(10):4592–9.