Image Anal Stereol 2002;21:13-18 Original Research Paper
SEMIAUTOMATIC DETECTION OF TUMORAL ZONE
Ezzeddine Zagrouba and Walid Barhoumi
Faculte des Sciences de Tunis, Groupe de Recherche Images et Formes de Tunisie 1060 le Belvedere Tunis, TUNISIA. e-mail: walid.barhoumi@enim.rnu.tn
(Accepted February 19, 2002)
Laboratoire CRISTAL,
ABSTRACT
This paper describes a robust method based on the cooperation of fuzzy classification and regions segmentation algorithms, in order to detect the tumoral zone in the brain Magnetic Resonance Imaging (MRI). On one hand, the classification in fuzzy sets is done by the Fuzzy C-Means algorithm (FCM), where a study of its different parameters and its complexity has been previously realised, which led us to improve it. On the other hand, the segmentation in regions is obtained by an hierarchical method through adaptive thresholding. Then, an operator expert selects a germ in the tumoral zone, and the class containing the sick zone is localised in return for the FCM algorithm. Finally, the superposition of the two partitions of the image will determine the sick zone. The originality of our approach is the parallel exploitation of different types of information in the image by the cooperation of two complementary approaches. This allows us to carry out a pertinent approach for the detection of sick zone in MRI images.
Keywords: FCM, fuzzy classification, MRI, segmentation, tumour brain.
INTRODUCTION
The objective of this work is to provide for the clinicians a semiautomatic system of diagnosis permitting the characterisation of the healthy and pathological matters from the digitalisation of the Magnetic Resonance Imaging (MRI) of the human brain. A great amount of work has focused on the treatment of the MRI brain images which is a complex task, considering the variability of the human brain and the complexity of the images. Among these works, we mention the approaches where the objective is the image segmentation (Joliot and Mazoyer, 1993; Brummer et al, 1993) and the approaches aiming at the classification of the image in matters (Tsao et al., 1992; Alyward et al, 1994; Vinitski et al, 1995). Our approach is mainly characterised by the integration
of the two methods realising complementary tasks. Indeed, we cooperate a system of segmentation in regions and a system of classification in fuzzy sets. First, the toboggan algorithm (Faiefield, 1990) is applied to the original image, in order to eliminate the noise often present in the medical images. Then, the image is split into a set of classes grouped in three matters (Grey, White and Other) and simultaneously segmented in regions. The class containing the sick zone is usually formed by several related components (regions). Next, an operator expert determines a germ (seed pixel) in the center of the sick zone as well as the size of the zone of interest (a window x x x). Finally, the extraction of the regions belonging to the intersection of the zone of interest and the tumour class allows to detect precisely the sick zone (Fig. 1).
Classification
Pre-treatments
Superposition —I
Segmentation
Operator Expert
Fig. 1. Proposed system for the sick zone detection.
13
Barhoumi W et al: Detection of tumoral zones
METHODS
FUZZY CLASSIFICATION
The classification in Fuzzy C-Means (FCM) is a generalisation of the classical clustering to the fuzzy sets domain. The advantage of fuzzy sets is the better modelling of uncertainty and ambiguity. This characteristic of fuzzy sets allows to postpone the classification of a pixel until more information can be used to make the final classification. Within the segmentation context, feature vectors are defined relatively to each pixel in the image. Using these feature vectors, FCM will create a fuzzy membership partition that gives a fuzzy membership for each pixel vector within each class. In fact, given a priori knowledge of the classes number c, FCM assigns a c-dimensional fuzzy membership vector to each feature vector and thus to each pixel. Then, FCM will cluster patterns so as to minimize the quadratic error between the fuzzy class centres and the pixels. Thus, FCM minimizes the energy function Jm(U, V;X) givenX to resolve the optimisation problem (n) (eq. 1).
create a random matrix Ut=0).
{n1 x n2 x c)  membership
o    Step 1: calculate the matrix Vt of the class centres
suchas:Vt=F(Ut) = fk(Ut))1<k<c
o    Step 2: calculate Ut+1) such as: Ut+1) = G(Vt) =
1<i<n11<j<n2,1<k<c
o    Step 3: if ||Ut+1)-UtA < e then terminate, else t = t+1 and go back to Step 1.
In  our  case,  the  functions fk  and gkj  are  given respectively by the equations (eq. 2) and (eq. 3):
n 1 n 2(uitjkmX i
ij
(2)
IKuitj(k)m
i1j=1
ij
{*):{
m i n 1 n 2 c((uijkm(
X ij-v k
A):
Vi,Vj,k:    uitj1' (k)=gki j(Vt>)
(2/m-1)
C = (Ue (0,1)n 1 n 2c c uij(k) = 1,Vi,Vj,
i
¦q1-
X  .-vt ij     k
Card
X i j-v q
ifI it> = <b
hj     r
otherwise
and    0<rni 1X^u  jk)<N,Vk).                         f
j                     i          withIitj = k/1<k
(3)
<c,
with N: the number of pixels to classify in the image Im(n1;n2) (N = n1 x n2),
p: the number of the chosen attributes to characterise a pixel,
X ij-vk
t
0J,^I is the
the
Xij G Rp where 1 < i < n1, 1 < j < vector representing the pixel Im(i,j),
vk G Rp where 1 < k < c: the centroid vector of the class k,
uij(k) G [0,1]: the degree of membership of the vector Xi j in the fuzzy class k,
m G ]1,+°°[: a weighting exponent controlling the fuzziness degree of the class,
U = (ui(k)) where 1 < i < n1, 1 < j < n2,1 < k <c: the three dimensional fuzzy matrix partition
\\X< norm A.
ij
¦v.
the distance between Xi   and vk in
Steps of the FCM algorithm according to Bezdek (1981) are the following:
o    Initialisations:
-   set parameters: c, m, A and the termination constant £;
characteristic function of the set I and CardI) is the cardinal of I.
Next, we study the FCM parameters in order to improve the results and to reduce the calculating time.
Feature extraction: several characteristics relative to the pixels can be used as vector components. Among them, we mention the intensity, the spatial location, the average and the variance of a w x w window centred in the considered pixel. The clustering goal is to produce a feature space quantifying the similarity between pixels belonging to the same tissue type. If we consider an homogeneous area, the variance calculation gives higher values for pixels near the area boundary than those in the area interior. Consequently, the considered homogeneous area will be split into two classes. Thus, the variance is not a good component. Moreover, our experimentations show that the use of the average and the spatial location does not improve the results (Fig. 2 b, c, d) and increases the calculating time. Therefore, only the intensity (p = 1) is used as a component of feature vector Xi j. Then, the classification of a pixel depends only on its intensity [Xi j = GreyLevel(Im(i,j))] what allows to reduce the data set.
14
Image Anal Stereol 2002;21:13-18
Class number: since the number of clusters is an input to the clustering algorithm, it needs to be determined before hand. Several methods have been suggested in the literature for determining dynamically an optimal value of c (Gath and Geva, 1989; Yager and Filev, 1994). separate the MRI brain scans into Grey, White and Other matter, then the number of classes must be 3 at a minimum. Due to the intensity distributions of the 3 matters in the human brain scans, the minimum number of classes is practically 4. The reason for this minimum is that the Cerebral Spinal Fluid (CSF) matter generally has a higher intensity than white brain matter for MRI images. It is also possible to perform fuzzy clustering for higher value of c and interpret more than one of these classes as Grey, White or Other matter. In the case of c = 4, the classes have been ordered according to the value of their centroids. Letting c1 signify the class whose centroid has the lowest value and c4 signify the class with the highest centroid value, the results were interpreted as follows:
-    c1 and c4 were classified as Other matter
-    c2 was classified as Grey matter
-    c3 was classified as White matter.
Fig. 2. Classifications with various p and c values (m = 3.5): (a) original image, (b) p = 1 {intensity} c = 4, (c) p = 2 {intensity, variance] c = 4, (d) p = 3 {intensity, variance, space location} c = 4, (e) p = 1 {intensity}, c = 5, (f) p = 1 {intensity} c = 6.
Nevertheless, it is difficult to say that there does not exist an even higher value of c for which an interpretation can give better results than when c = 4 (Fig. 2 b, e, f). However, increasing c has the added expense of increasing the memory requirements for the algorithm.
Norm: the norm matrix A determines which type of distance is used, so that matrix A must be positive definite. The three most common distance norms are the Euclidean, the Diagonal and the Mahalanobis. The latter being frequently used because of its consideration of data distribution although it increments the FCM cost. For this, (Idrissa and Acheroy, 2000) apply Kahunen-Loeve transformation in order to reduce the algorithm cost.
Coefficient m: m is a weighting exponent that determines the degree of fuzziness of the fuzzy clusters. Our experimentations show that a high value of m makes the sets boundaries ambiguous (Fig. 3). That is to say, a high value of m delays the classification task.
Fig.  3.  Classifications with various m values:  (a) original image, (b) m=3.5, (c) m=6.5, (d) m=8.5.
In Fig. 4 a, we applied the FCM algorithm with various values of m, and for every value of m we determine the classification rate CR defined by:
CR = Card(8) ,    with
8 = I Im(ij)/ max ui (k) > 1 - eCR \ .   (4)
We conclude that the bigger m is, the smaller CR is, and therefore the decision of classification is still doubtful. However, a great value of m decreases the calculating time T (Fig. 4 b).
as.-------------------------,------------------	
-        \	<<«
•1	
2            i            6             9           10           12           H	
; 1		
\ <w	(0	
	\/	
2            4	B            10	n        h
Fig.   4.   Variation   according   to   m   of:   (a)   the classification rate (CR), (b) the calculation time (T).
15
Barhoumi W et al: Detection of tumoral zones
In order to take into account the previous deductions concerning the variation of CR and T relatively to m, we choose a great value of m and reduce it during the algorithm according to the classification rate. At every iteration t, the more CR increases, the more we decrease m (eq. 5). In fact, this adaptation allows to fairly reduce the calculating time of 10.5% while making better the quality of clustering (i.e.CR~N).
if   CR > q 2t
then    m := m — q
3 ¦
(5)
In practice, we empirically define: q 2 t = (2 + t)N/10 and q3 = 0.25.
As the random initialization of the partition matrix U 0) influences the outcome of the classification algorithm, we choose to initialise the class centroid matrix V(0) (eq. 6) so that the centroids evenly span the full range of input data and then the algorithm begins from the step 2. This is done in order to guarantee a better distribution of the class centres.
Vhe{1,2,-,p},Vke{1,2,-,c}:
vk(h)=kE(1       max      Xjh) -
\c  1<i<n1,1<j<n2      'J
1<<min<< X ij h])- (6)
Once the FCM algorithm is finished and the final matrix partition Uf is determined, the final classification will be hard. In fact, each pixel Im(i,j) will be clustered in the class Ck (k* e {1,--- ,c}) according to:
Af)
fj(k*
max uVlik)
1<k<c  i<j v  !
(7)
HIERARCHICAL SEGMENTATION BY ADAPTIVE THRESHOLDING
This stage consists in providing a segmentation of the image into regions. We use our algorithm of hierarchical segmentation by adaptive thresholding (Barhoumi et al., 2000) which provides better results than the classical techniques such as Split-and-Merge (Zagrouba et al., 1994), Region Growing (Haralick and Shapiro, 1985) and Thresholding (Zagrouba et al., 1998). Our approach consists in using several image representations with levels of increasingly fine details (coarse-to-fine). Thus, we build a tree of encased segmentations.
Initially, we consider the image as a whole region to be the hierarchy root. At each level, each node (region) of the tree is either a leaf (indivisible region)
or it has n children corresponding to its splitting into sub-regions {R1,R2,--- ,Rn}. This process continues until no further splitting occurs. The final segmentation result of this top-down method is obtained by the union of the indivisible regions (terminal nodes in the tree structure) found at every level.
For this, we consider that a region is defined as the set of interior points delimited by closed contours characterised by a higher gradient value (Canny, 1986). Then, we classify as interior points to a region, the points of lower gradient than a threshold S*, and the other points will delimit the region. Thus, the splitting of a region R amounts to find the suitable threshold S*. Then, indeed S* is determined, the partition P of R (P = {R1R2,--- ,Rn}) will be the set of related components having lower gradient than S*.
Automatic determination of the threshold: as segmenting a region R consists in dividing it into sub-regions that we hope more homogeneous, we evaluate the partition P = {R1,R2,--- ,Rn} of R with a monotonous function FR (eq. 8) based on the calculation of variance which is a good feature characterising the regions homogeneousness.
Fr(s) = Varianc RiR) å Size(Ri) Variance Ri  _   (8)
Therefore, for every region belonging to this tree and appearing in previous steps with a threshold SR, the segmentation threshold S* will be the biggest in the interval [0,SR] for which the leap of the function FR will be upon a threshold f (eq. 9). This corresponds to the appearance of at least one new significant sub-region.
S* = Sup{Se0,SR]/FR(S) >f.
(9)
Then, it remains to decide when a branch of the hierarchy must stop ramify what corresponds to the inability of the correspondent region to produce significant sub-regions.
Termination criterion: if the distribution of the gradient values on a region R is almost-uniform, then R is considered as indivisible. To estimate the distribution of the gradient, we use the gradient histogram hRgrad of R. Then, we calculate relatively to R the entropy ER which provides an idea on the partition of the histogram values along the thresholds axis (eq. 10).
ER = -åhgRrad(s)log(hgRrad(s)),         (10)
s0
16
Image Anal Stereol 2002;21:13-18
with   hRrad{S) = nR S/å SR 1nRs)   and   nR(s) =     the sick zone. This class CK is concentrated and limited
in a unique zone of the brain (CK C V(ger0)). However, the superposition of CK with the regions map allows to detect precisely the three regions composing the sick zone. Then, Fig. 6 shows the results of the application of this approach in the case of an another type of tumour. In this case, the class Cv is not concentrated in
Card N (i, j) G R /  grad Im (i, j)   = s > j.
Knowing that ER is minimal when the frequencies are nil almost everywhere except for a value and is maximal when the frequencies distribution is uniform, this notion of entropy allows to decide whether a region is divisible. Indeed, if ER is maximal ([ER > log(S)]) then R is indivisible and it is considered as a leaf of the hierarchy.
SUPERPOSITION OF THE TWO PARTITIONS
Our experimentations showed that the tumoral zone is represented by a set of regions belonging to the same class of matter CK which is not entirely sick. In fact, some healthy pixels have a similar grey level to sick pixels. So, this first approximation of the tumour (class CK) need to be refined. Given the class map generated by the FCM algorithm, the selection of a seed pixel ger0 in the tumoral zone of the image by an operator expert determines the class CK. Moreover, the operator expert defines a zone of interest V(ger0) (a window x x x), centred on ger0, in order to limit the searching area of the sick zone. Moreover, healthy pixels still remain in the zone of search (CK n V(ger0)). But, the gradient information allows the segmentation algorithm to isolate sick pixels from their neighbouring healthy ones. Given the regions map generated by the hierarchical segmentation, the regions belonging to the intersection of CK and V(ger0) form the pathological zone Zmal (Figs. 5 and 6).
Formally, if FS is the function of segmentation
r
verifying  F S(Im) = {R1 R 2,-,Rr}   with   [JRi C
i=1
Im and FC the function of classification verifying
c
F C(Im) = {C15C2,--- ,Cc}   with   [JCi = Im, then it
i=1
exists an integer K G {1,--- ,c} verifying ger0 G CK and consequently Zmal C CK. Finally, the sick zone will be Zmal = [JRt such as: Rt C (FS(Im) nCKnV(ger0)).
RESULTS
Fig. 5 shows the successive steps of the application of our global approach on an MRI tumoral brain (p=1,c = 4, m 0 = 6.5, f = 0.01, eCR = 0.2). The selection of the germ ger0 and the application of the FCM algorithm allows to determine the tumour class CK which can be considered as a first approximation of
k -a limited area of the brain (CK not included in V(ger0)).
In fact, the class CK is formed by the tumoral zone
and a set of pixels distributed along the brain boundary
(Fig. 6 b). Then, the zone of interest is necessary to
put out these boundary pixels of the searching area.
Finally, the superposition step allows the detection of
the sick zone as an aggregate of small regions.
V (ger,)
Fig. 5. (a) original image, (b) classified image, (c) segmented image, (d) superposition.
V (gM-o)
Fig. 6. (a) original image, (b) classified image, (c) segmented image, (d) superposition.
CONCLUSION AND PERSPECTIVES
The cooperation of the classification and the segmentation algorithms allows to carry out a robust approach of the sick zone detection in the brain MRI images. In fact, our algorithm of classification is an improvement of the classical FCM algorithm allowing to get better results and reduce the calculating time. This classification could serve as a relatively quick precursor to the separation of anatomical structures. Moreover, our segmentation algorithm produces better results than the traditional methods. Then, let us note that our approach can be applied to other types of images with the aim of detecting the zones
17
Barhoumi W et al: Detection of tumoral zones
of interest. Finally, we can perform our approach by applying other combinatorial methods (Taboo, Monte Carlo, Genetic algorithms, Neural networks,...) to solve the problem p and by paralleling the segmentation algorithm considering the independence of the regions processing in the built tree which reduces the calculating time.
ACKNOWLEDGMENTS
We sincerely thank Mr. Yevgeniy Kuznetsov, Neurological Institute and Hospital, Montreal, Canada; who allowed us to recover the MRI images appearing in this article.
REFERENCES
Alyward S, Coggins J, Cizaldo T, Andereasen N (1994). The effects of magnetic resonance image inhomogeneities on automated tissue classification., TR94-021, University of North Carolina, Chapel Hill.
Barhoumi W, Zagrouba E, Ghorbel F (2000). Hierarchical region segmentation and matching for stereoscopic images, World Multiconference on Systemics, cybernetics and Informatics, Orlando, USA. V:16-20.
Bezdek JC (1981). Pattern recognition with fuzzy objective function Algorithms. New York: Plenum Press.
Canny J (1986). A computational approach to edge detection. IEEE Trans Pattern Anal, 8(6):679-98.
Faiefield J (1990). Toboggan contrast enhancement for contrast segmentation. 10th International Conference on Pattern Recognition. Atlantic City.
Gath I, Geva B (1989). Unsupervised optimal fuzzy clustering. IEEE Trans Pattern Anal, 2(7):773-81.
Haralick RM, Shapiro LG (1985). Image segmentation techniques. CVGIP-29-1, 100-32.
Idrissa M, Acheroy M (2000). An iterative fuzzy c-means algorithm-application to Landcover classification. World Multiconference on Systemics, cybernetics and Informatics, Orlando, USA. V:224-6.
Joliot M, Mazoyer BM (1993). 3D segmentation and interpolation of magnetic resonance brain images. IEEE Trans Med Imaging, 12(2):269-77.
Brummer ME, Mersereau RM, Einsner RL, Lewine RRJ (1993). Automatic detection of brain contours in MRI data sets. IEEE Trans Med Imaging, 12(2):153-65.
Tsao ECK, Lin WC, Chen CT, Bezdek JC, Pal NR (1992). A neural network system for medical image understanding. Proceedings of the 5th Florida Artificial Reserach Symposium, Florida AI Research Society.
Vinitski S, Gonzalez C, Burnett C, Seshagiri S, Mohamed F, Ortega H, Faro S (1995).Tissue segmentation in MRI as an informative indicator of disease activity in the brain. 8th ICIAP, Sanremo, Italy.
Yager RR, Filev DP (1994). Approximate clustering via the mountain method. IEEE Trans Syst Man Cyb, 24(8):1279-389.
Zagrouba E, Buvry M, Krey C (1994). A rule-based system for region segmentation improvement in stereovision. IS&T-SPIE International Technical Conference of San Jose, Vol 21-28, 357-67. USA.
Zagrouba E, Hedidar T, Jaoua A (1998). A robust and unified algorithm for indoor and Outdoor Scenes on region segmentation. The 12th IEA-AIE, Cairo, Egypt. 570-9.
18