
Research Article


Medical Image Segmentation Using Enhanced HoshenKopelman Algorithm 

K. Kies
and
N. Benamrane



ABSTRACT

This research present a method of 2D and 3D segmentation based on the
Enhanced HoshenKopelman algorithm and its extension to nonlattice structure.
The main feature of this method is to combine a merging strategy of a
region growing algorithm with the multiple labeling technique of the EHK
algorithm for regular and nonregular lattice. An efficient reconstruction
algorithm is then applied to the set of edge points for the obtained segmented
regions. The latest uses 3D Delaunay weighted triangulation. The combination
of these already known algorithms makes the proposed approach very fast,
efficient and appropriate to volumetric segmentation and particularly
the anatomical structures in MRI and CTScan images.







INTRODUCTION
Segmentation is an indispensable tool in medical image processing
and is useful in many applications such as early pathology detection,
precise tumor localization, simulation and planning of chirurgical interventions,
etc. The principal aim of segmentation is the extraction of significant
objects present in the image either by partitioning this image into connected
semantic regions or by extracting one or many specific objects from this
image.
Proposed approaches attempting to solve the segmentation problem are
numerous and extremely varied. The choice of a method depends on the nature
of the image, what primitives are to be extracted and what are the usage
constraints (algorithm complexity, real time requirement, available memory,
etc...). This variety of possible methods shows that the segmentation
problem is far from being solved and there exists no single technique
satisfying all the image types.
In the literature, two major classes of segmentation algorithms are present:
contour approach based on the estimation of the locale discontinuities
of the gray level function and region approach based on homogeneous region
extraction. Among the later, a fundamental algorithm is region growing
which has received considerable attention from 2D and 3D image analysts
and researchers. The principle of the technique is to define some criteria,
called homogeneity criteria, which allow regrouping connected image pixels
into larger regions. These criteria are often based on differential surface
characteristics of the 3D object and a stable segmentation in quadric
pieces is achieved (Monga and Bricault, 1994). They can also be obtained
from some evaluation function (Revol and Jourlin, 1994) or by RMSapproximation
of image points by simple polynomial surfaces such as planes or quadrics
(Djebaili, 1998). Region growing method has also been used for 3D medical
image reconstruction where the voxel gray level defines some density of
the volume to reconstruct as well as a regrouping criterion.
Another approach that is fundamentally geometrical which is mathematical
morphology mainly consists on comparing objects under analysis to other
objects with known shapes. Morphologic operators applied to 3D images
are described in Agnus et al. (2000) where temporal image sequences
considered as 3D images are segmented. Based on mathematical morphology,
threedimensional segmentation has been applied to anatomical structures
in MR images on large databases (Bueno et al., 2000).
In the last two decades, a promising mathematical framework based on
variational models and partial differential equations has been developed
and used to solve the segmentation problem. There exist two major categories
of such models: parametric models introduced in Kass et al. (1988)
and geometric models based on curve evolution theory and level sets (Caselles
et al., 1993; Malladi et al., 1995). They allow segmentation
and object reconstruction by providing compact geometric representations
where evolution laws similar to physical mechanics laws are introduced
such as forces, physical constraints and friction (Cohen, 1992; Lachaud,
1998; Montagnat, 1999; Montagnat et al., 2001). Procedures using
deformable models have been widely used for segmentation since they are
well suited for object extraction from temporal and volumetric data from image sequences such
as IRM, CT scanned, X Tomography, ARM and Echography (Cohen, 1992; Montagnat,
1999; Delingette and Montagnat, 1999). These procedures have been made
more efficient by using 3D pyramid images (Lachaud and Montanvert, 1999)
and such hierarchical structures allow coarse to fine type segmentation.
Region and edge information has also been integrated into these models
(Pardo et al., 2001).
In some applications, the segmentation is viewed as a pattern recognition
problem. Thus, the task of delimiting regions is associated to a high
level ability to classify regions in particular semantic classes. These
approaches are widely used for segmenting cerebral MRI; a review is given
in Bezdek et al. (1993). Among these, we can cite: statistical
methods using Markov fields and Bayesian classification (Jaggi, 1998),
the non parametric methods such KNN algorithms introduced by Duda and
Hart (1973), Parzen windows, fuzzy logic based methods (Bezdek et al.,
1997) and the so called Biologists methods: neural networks, genetic algorithms,
etc. (Takatsuka and Jaivis, 2001; Kobashi et al., 2001; Kim et
al., 2001).
In quest of better efficiency and precision various hybridations of these different
techniques have been successfully attempted (Kapur et al., 1996; Geraud,
1998; BenoitCatlin et al., 2000; Helliera and Barillot, 2001; Rey et
al., 2002).
In this research, a method for 2D and 3D image segmentation based on
the HoshenKopelman algorithm is proposed. The main feature of this method
is to combine a merging strategy of a region growing algorithm with the
multiple labeling technique of the EHK algorithm for regular and nonregular
lattice. This procedure is very efficient and fast and may be adapted
to the segmentation of anatomical structures in MRI and CTS can images
which are in general quite large. An efficient reconstruction algorithm
is then applied to the set of edge points for the obtained segmented regions.
The latest uses 3D Delaunay weighted triangulation. The later is very
efficient and fast approach makes it possible to adapt it to the segmentation
of anatomical structures in MRI and CT Scan images.
THE HOSHENKOPELMAN ALGORITHM
Connected component analysis of grey level image divides this image
into regions by labeling uniquely each pixel to a defined region. The
regions discriminated by the connected component depend on the definition
of neighborhood. Several definitions of neighborhood exist, however the
most commonly used is the 4connected. Basically, the algorithms for connected
components assign new labels to the first pixel of each component and
attempt to propagate these labels through the neighbors. In image processing,
Rosenfeld and Pfalz (1966) described an algorithm based on the classical
component algorithm for graphs which is known as the classical algorithm
for connected component analysis. The algorithm makes two passes through
the image. The first pass assigns and propagates labels, while at the
same time it creates an equivalence class table. The second pass reassigns
the labels that belong to the same class using the equivalence class table,
such that all equivalence is assigned with a unique label. The final image
has its connected components indexed with a unique label. While being
relatively fast this algorithm requires a large amount of memory to store
the equivalent class table, which depends on the image size and on the
complexity of the images. The most recent algorithms, for connected component
analysis, attempt to overcome the fundamental problems of either being
too slow or requiring a large amount of computer memory.
In physics connected component analysis is known as cluster analysis
and it is widely used for problems related to the percolation phenomena.
In mid 1970s, Hoshen and Kopelman (1976) published an ingenious and efficient
cluster analysis algorithm, which is known as the HoshenKopelman (HK)
algorithm. The time complexity of this algorithm is linear and requires
small computer memory size. It remains an important breakthrough in statistical
physics, since it allows the study of cluster size and the percolation
problem through computer simulations. Recently, Hoshen (1998) has developed
the Enhanced HoshenKopelman (EHK) algorithm, which is an improvement
toward cluster analysis. It is capable to yield various shape parameters
and still preserves the computational time and memory space complexities
of the original HK algorithm. The key idea for the efficiency of the HK
algorithm, which is also essential in the EHK algorithm, lies on the concept
of multiple labeling technique. The multiple labeling technique assigns
different labels to pixels belonging to the same cluster. Unlike the classical
connected component algorithm, which stores an equivalence class table
and the equivalence are executed by a second pass over the image. The
HK algorithm generates the proper label by assigning negative indexes
to equivalent labels and recursively searching for the final label.
In the statistical physics terminology, an image is referred to as a lattice,
the connected components are called clusters and a pixel is known as a site.
The algorithm begins by scanning the lattice, pixel by pixel or site by site,
from up to down and right to left direction (and back to front for 3D images).
SEGMENTATION
Grey level segmentation may be achieved using standard region growing
technique which is a merging process at various levels. Usually, two or
more steps are used; the first step (step 1) operates at the pixel level
which is referred as a presegmentation step. At this stage the lattice
is regular. The further segmentation merging steps consist on whether
regions obtained from step 1 may be merged or not. The first two segmentation
steps may make use of uniformity predicates associated to the characteristic
attributes of region say Ri are as follow:
• 
P_{1}(Ri) = [(MAXi MINi) < Threshold_{1}] where
MAXi and MINi are respectively the minimum and maximum gray level
in region Ri 
• 
P_{2}(Ri) = [V(Ri) < Threshold_{2} ] where V(Ri)
is the gray level variance of points in region Ri. We have Threshold_{2}
= (Threshold_{1})^{2} 
The thresholds may be chosen in any standard manner, this issue is not
discussed here. Other merging steps may make use of other predicates such
as shape or texture.
In this research, we propose using a variation of the EHK algorithm for
grey level image segmentation. While the main features of the algorithm
are preserved, the notion of non occupied site does no longer apply here
since all sites (or pixels) are labeled. The background is also a region.
For step 1, the standard HK algorithm is used since at the pixel level,
the lattice is regular. The clusters obtained correspond to the regions
of the presegmentation step. After the completion of step 1, a region
adjacency graph may be constructed. This graph corresponds to a non regular
lattice. This latter may become quite large for very large images. In
order to perform step 2 and further merging steps in an efficient manner,
we propose and make use of EHK algorithm to nonlattice environments as
described by AlFutaisi and Patzek (2003).
Relative to other region growing segmentation methods which are at best
O(nlogn) in time complexity (n being, the number of pixels or voxels),
this method is linear in O(n). In Fig. 1, we have plotted
the computational time in seconds for random binary 3D images of sizes
ranging from 100^{3} to 300^{3} (a 2.6 GHz Intel Pentium
IV with 2 GB DDR400 RAM PC was used).
Two more features make this method very attractive.

Fig. 1: 
Cpu times as function of 3D image size in voxels; our
method (+), standard region growing method (x) 
Firstly, the memory
requirement depends solely on the maximum number of labels needed and not on the size of the image. Secondly,
this method can be highly parallelized as shown in Tiggemann`s work on
percolation (Tiggemann, 2001) where system sizes as high as 20000^{3}
where considered.
RECONSTRUCTION
Once done with previous steps, each region of the segmented 3D image
is delimited with a set of edge points from which we can get the region
surface using appropriate reconstruction method. In practice, even for
modest regions, the number of points can be quite large (for example a
20 pixel diameter cylinder across 10 slices counts some 600 points) while
large regions may require more than 100,000 points. Following (Kies et
al., 2005), were efficient 3D Delaunay triangulation (3DT) algorithm
is used. The outcome of 3DT is a set of pair wise linked edge points,
in which each linked triplet forms a triangular face and each face quadruplet
forms a tetrahedron. The outer surface or the crust of the region of interest
is a subset of the set of all triangular faces resulting from 3DT. If
the crust is convex that the solution is straightforward: it is the Convex
Hull. To tackle the problem of the Concave Hull, much elaborated algorithms
have been developed.
In the reconstruction process, we used principally Nina et al.
(2001). We have also used Patrice Koehl`s (2002), TETRAFOR which has very
efficient 3D Delaunay Triangulation (Weighted). Once the crust faces obtained,
there is one last step which consists on correctly orienting the normal
of each face toward the outside of the object. The success of this last
step allows for a 3D rendering with no defects during final visualization
of the results.
RESULTS
First of all and in order to validate our results, we have tested
this method on various test 3D images containing various simple objects
(cubes, spheres and cylinders). We have also tested random binary images
of various sizes and various ratio black/white. We then segmented a 256
slices 256x256 scaled to 8 bits MRI Brain 3D image and extracted and reconstructed
the brain surface as the largest region. Figure 2 shows
several slices.

Fig. 2: 
A set of MRI Brain Image slices 

Fig. 3: 
A bottom view of the whole brain reconstructed from
50% of all the edge points 
Figure 3 represents the outer surface reconstructed using a subset of edge points for the largest segmented
region representing the whole brain.
Similar results were obtained from 256x256x256 CT Scanner image of the
skull of a Monkey. Figure 4 shows a sample set of slices
used.
The 3D rendering after using half of all edge points as input for powercrust
is shown in Fig. 5.

Fig. 4: 
A set of CT scanner image slices 

Fig. 5: 
3D rendering of crust obtained for the region representing
the Skull of a Monkey 
CONCLUSIONS
We have presented a new and very efficient method of 3D segmentation
base on both the Enhanced HoshenKopelman algorithm for both lattice and
nonlattice environment and applied it to some test images. We have illustrated
our results by reconstructing the largest region using POWERCRUST. The
efficiency of the proposed segmentation approach makes it possible to
tackle very large medical images with very small CPU time. An ongoing
work is underway for further benchmarking. Also, an implementation of
this method on parallel machines is in preparation where better performance
is expected.

REFERENCES 
1: Agnus, A., C. Ronse and F. Heitz, 2000. Segmentation spatio temporelle morphologique de saquences d'images. LSIIT, UPRESA 7005, University de Strasbourg.
2: AlFutaisi, A. and T.W. Patzek, 2003. Extension of the hoshenkopelman algorithm to a nonlattice environment. Physica A, 321: 665678. Direct Link 
3: BenoitCattin H., T. Zouagi and C. Odet, 2000. Une vision fonctionnelle de la segmentation d’images. CREATIS, Unité de Recherches CNRS (UMR5515).
4: Bezdek, J.C., L.O. Hall and L.P. Clarke, 1993. Review of MR image segmentation techniques using pattern recognition. Med. Phys., 20: 10331047. Direct Link 
5: Bezdek, J.C., L.O. Hall, M.C. Clark, D.B. Goldgof and L.P. Clarke, 1997. Medical image analysis with fuzzy models. Stat. Methods Med. Res., 6: 191214. PubMed  Direct Link 
6: Bueno, G., O. Musse, F. Heitz and J.P. Armspash, 2000. Threedimensional analysis of anatomical structures in MR images on large databases. University Louis Pasteur (Strasbourg I), LSIIT. CNRSUPRES A7005.
7: Caselles, V., F. Catt, T. Coll and F. Dibos, 1993. A geometric model for active contours in image processing. Numerische Math., 66: 131. CrossRef 
8: Cohen, I., 1992. Modèles déformables 2D et 3D: Application à la segmentation d’image médicales. Ph.D Thesis. Université de Paris IX Dauphine.
9: Delingette, H. and J. Montagnat, 1999. Extension des maillages simplexes pour la segmentation images 3D et 4D. Projet Epidaure.
10: Djebaili, M., 1998. Segmentation des Images de profondeurs utilisant les analyses multiresolutions. Thyse LIGIM.
11: Duda, R.O. and P.E. Hart, 1973. Pattern Classification and Scene Analysis. 2nd Edn., John Wiley and Sons, New York, ISBN: 0471223611, pp: 482.
12: Geraud, T., 1998. Segmentation automatique des structures cyybrales internes en imagerie par rysonance magnytique tridimensionnelle. Ph.D. Thesis, University de laENST Paris.
13: Helliera, B.P. and C. Barillot, 2001. Segmentation of brain 3D images using level sets and dense registration. Med. Image Anal., 5: 185194. Direct Link 
14: Hoshen, J. and R. Kopelman, 1976. Percolation and cluster distribution. I. Cluster multiple labeling technique and critical concentration algorithm. Phys. Rev. B, 14: 34383445. CrossRef  Direct Link 
15: Hoshen, J., 1998. On the application of the enhanced HoshenKopelman algorithm for Image analysis. Pattern Recog. Lett., 19: 575584. CrossRef 
16: Jaggi, C., 1998. Segmentation par méthode markovienne de l’encéphale humain en imagerie par resonance magnétique: Théorie, mise en oeuvre et evaluation. Ph.D Thesis. Université de Caen, France.
17: Kapur, T., W.E.L. Grimson, W.M. Wells and R. Kikinis, 1996. Segmentation of brain tissue from magnetic resonance images. Med. Image Anal., 1: 109127. Direct Link 
18: Kass, M., A. Witkin and D. Terzopoulos, 1988. Snakes: Active contour models. Int. J. Comput. Vision, 1: 321331. CrossRef  Direct Link 
19: Kies, K., N. Benamrane and A. Benyettou, 2005. 3D medical image segmentation and surface modeling using the power Crust. Inform. Technol. J., 4: 377381. CrossRef  Direct Link 
20: Kim, E.Y., S.W. Hwang, S.H. Park and H.J. Kim, 2001. Spatiotemporal segmentation using genetic algorithms. Pattern Recog., 34: 20632066. CrossRef  Direct Link 
21: Kobashi, S., N. Kamiura, Y. Hata and F. Miyawaki, 2001. Volumequantizationbased neural network approach to 3D MR angiography image segmentation. Image Vision Comput., 19: 185193. CrossRef  Direct Link 
22: Koehl, P., 2002. A program for computing 3D weighted Delaunay triangulations (full source code under LGPL). http://csb.stanford.edu/~koehl/ProShape/tetrafor.tar.gz.
23: Lachaud, J.O., 1998. Extraction de surfaces à partir d'images tridimensionnelles: Approche discrète et approche par modèle déformable. Ph.D Thesis. Université Joseph Fourier, Grenoble, France.
24: Lachaud, J.O. and A. Montanvert, 1999. Deformable meshes with automated topology changes for coarse to fine dimensional surface extraction. Med. Image Anal., 3: 187207. CrossRef 
25: Malladi, R., J.A. Sethian and B.C. Vemuri, 1995. Shape modeling with front propagation: A level set approach. IEEE Trans. Pattern Anal. Mach. Intell., 17: 158175. CrossRef  Direct Link 
26: Monga, O. and I. Bricault, 1994. From volume medical images to quadric surface patches. R.R. No. 2380. INRIA.
27: Montagnat, J., 1999. Modèles déformables pour la segmentation et modélisation d’images médicales 3D et 4D. Ph.D Thesis. Nice SophiaAntipolis.
28: Montagnat, J., H. Delingette and N. Ayache, 2001. A review of deformable surfaces: Topology, geometry and deformation. Image Vis. Comput., 19: 10231040. CrossRef 
29: Nina, A., S. Choi and R. Kolluri, 2001. The power crust. Proceedings of the 6th ACM Symposium on Solid Modeling, June 48, 2001, Ann Arbor, Michigan, United States, pp: 249260.
30: Pardo, X.M., M.J. Carreira, A. Mosquera and D. Cabello, 2001. A snake for CT image segmentation integrating region and edge deformation. Image Vision Comput., 19: 461475. CrossRef 
31: Revol, C. and M. Jourlin, 1997. A new minimum variance region growing algorithm for image segmentation. Pattern Recog. Lett., 18: 249258. CrossRef 
32: Rey, D., G. Subsol, H. Delingette and N. Ayache, 2002. Automatic detection and segmentation of evolving processes in 3D images: Application to multiple sclerosis. Med. Image Anal., 6: 163179. Direct Link 
33: Rosenfeld, A. and J.L. Pfalz, 1966. Sequential operations in digital picture processing. J. ACM., 13: 471494. CrossRef 
34: Takatsuka, M. and R.A. Jaivis, 2001. Encoding 3D structural information using multiple selforganizing feature maps. Image Vision, 19: 99118. Direct Link 
35: Tiggemann, D., 2001. Simulation of percolation on parallel computers. Ph.D Thesis, University of Cologne.



