Robust c-prototypes algorithms for color image segmentation

In this paper, we present a modified clustering algorithm to segment color images. The proposed technique is based on the robust c-prototypes algorithm with some modifications in its objective function and updating equations. The JK chromatic subspace of the IJK color space is used to segment color images. The algorithm performance is tested on real images with natural artifacts that make the segmentation process difficult. Additionally, a remote sensing image is segmented to demonstrate that the proposed algorithm can be used in real applications. Simulation results indicate that the proposed method is more effective than others proposed in the literature in terms of objective and subjective criteria.

Clustering is the process to discover and organize objects into subsets of similar objects called clusters. Pixel clustering in three-dimensional (3D) color space on the basis of color similarity is one of the current approaches in the field of color image segmentation. Colors in an image create dense clusters in the color space in a natural way and look like three pixel clouds representing clusters in the color space. Clustering is often seen as an unsupervised classification of pixels, and usually, a priori knowledge of the image is not used during the clustering process [1].

Different clustering techniques have been proposed to color image segmentation [2]. One of the most popular clustering techniques is the fuzzy c-means (FCM) technique [3]. In general, the FCM algorithm is a highly effective methodology to segment noise-free images, but in the presence of natural artifacts (noise, intensity, or color inhomogeneity in the regions, the regions with similar textures, shadows, object reflections, etc.), the FCM has two shortcomings that make it very sensitive: in its conventional nomenclature, it does not consider any spatial information in the image context [4–8], and the second one is that the objective function can be seen as a formulation of the least squares method, in which one tries to minimize the error between the feature vector and the vector with the centers of the groups. Outliers have a great effect during minimization since there is a quadratic function in the objective function of the FCM algorithm; so, it is necessary to use a quadratic function with the property of being less increasing, and thereby control the influence of atypical information.

Many researchers have taken FCM algorithm as the starting point and modify it to solve the first shortcoming. In [4], a spatial function is introduced into the membership function; this spatial function is the summation of the membership function in the neighborhood of each pixel under consideration. Ahmed et al. [5] modified the objective function of the FCM by adding a regularizer which considers the labels of all the neighboring voxels in an image to compensate the bias field effect of the data. Other authors such as Zhang and Chen [6] follow the same work of Ahmed et al. but their algorithms introduced the median- and mean- filtered images, which can be computed in advance and hence, can reduce the computation time. Zhao et al. [7] used an adaptive spatial parameter for each pixel that was designed to make the non-local spatial information of each pixel playing a different role to guide the noisy image segmentation. Liu and Pham [8] modified the FCM by adding a spatial penalty term into the objective function. The penalty term considers the fuzzy memberships from the neighbors of the current samples, often a pixel or voxel in 2D or 3D images. By this means, the clustering algorithm can take both the image-intensity information and spatial information into account and hence is more robust to noise and outliers producing better segmentation results. Other researchers have modified the FCM algorithm to get more robust segmentation hybrid algorithms by using other important approaches. For example, in [9], a method that combines the benefits stemming from the sound spatial coherency modeling capabilities of the Hidden Markov random field (HMRF) model and the enhanced flexibility obtained by the FCM algorithm was proposed. Siang and Ashidi [10] proposed an approach that analyzes the uniformity of the regions in color images by the histogram thresholding technique, and then the FCM is used to improve the compactness of the clusters forming uniform regions. A proposal based in the use of hyperplanes in the FCM algorithm and spatial constraints is presented in [11]; in this algorithm, a spatial regularizer is added into the fuzzy hyperplane-based objective function, but taking into account information of inherently spatial data. In [12], the authors presented an approach that relies on cluster center initialization and color quantization allowing faster and more accurate convergence such that it is suitable to segment very large color images. An analysis of different methods of initialization for the k-means algorithm for segmenting images is presented in [13]. Dae-Won et al. proposed in [14] an algorithm that extracts the most vivid and distinguishable colors referred as the dominant colors, the color points closest to these dominant colors are selected as the initial centroids in the FCM calculations using a fuzzy membership model between a color point and a reference color. Finally, an approach designed to promote an optimal initialization scheme of the FCM algorithm using a single point iterative weighted based upon the prior information was presented in [15].

A method for solving the second of the shortcomings of the FCM algorithm is based on the exploitation of the robust estimators especially the M-estimator [16]. This proposal was studied by Frigui and Krishnapuram [17], who designed their own robust estimator, based on special loss and weight functions. In this algorithm, the robust estimators allow the control of the influence of outliers in calculating the centers of the all groups. This paper is focused on the second shortcoming of the FCM algorithm. For this, we propose to use classical M-estimators and as the robust framework [18], states to use their loss and influence functions instead. We extent the proposed algorithm in a bidimensional clustering approach using the chromatic subspace in the IJK color space to segment color images [19].

2.1 Fuzzy c-means clustering

Fuzzy cluster analysis allows gradual membership of data points to clusters measured as degrees in [0, 1] providing the flexibility needed to express that data points can belong to more than one cluster. Additionally, the membership degrees offer a finer degree of detail of the data model and can express how a data point should belong to a cluster or not [16, 20].

Let X = {x ₁, …, x _n } be the set of given featured data and let c the number of clusters (1 < c < n) represented by the fuzzy set C _j (j = 1, …, c). Then, U _f  = (u _ij ) is called the fuzzy cluster partition of X if ${\displaystyle {\sum }_{i=1}^n{u}_{\mathit{ij}}<0}$, ∀j ∊ {1, …, c} and ${\displaystyle {\sum }_{j=1}^c{u}_{\mathit{ij}}=1}$, ∀i ∊ {1, …, n} hold. A fuzzy cluster model of a given data set X into c clusters is defined to be optimal when it minimizes the following objective function under the above constraints,

where ‖x _i  − c _j ‖² is the square of the Euclidean distance from feature vector x _i to the center of the class c _j and the parameter m ≥ 1 is a weighting exponent called the fuzzifier. The value of m determines the ‘fuzziness’ of the classification; when m = 1, there is a hard assignment of the data in the classes (k-means algorithm), on the other hand when m = 2 the data are gradually assigned with an interval of [0,1] to the different groups (FCM algorithm). The objective function J _f is alternately optimized using the membership degrees u _ij and the cluster centers c _j by setting the derivative of J _f with respect to the parameters equal to zero taking into account the constraint stated above. Finally, the equations for the two iterative steps that form the FCM algorithm are given as follows:

2.2 Robust m-estimator

M-estimators are a generalization of maximum likelihood estimation (MLE) and were proposed by Peter Hubert [18, 21–23]. Their definition is given by a robust loss function ρ(x) = ln (f(x)), connected with the probability density function f(x) for the sample data x _i , i = 1, …, n. The objective of M-estimators is to find an estimation $\widehat{\theta }$ of such that

The estimation of the localization parameter θ can be found by calculating the partial derivative of ρ with respect to θ introducing the influence function ψ(x, θ) = ∂ ρ(x, θ)/ ∂ θ,

Definition (5) is not fully equivalent to (4). However, it is used to find a solution for the minimization (4). Both definitions of the M-estimator are implicit, so are necessary iterative techniques for the calculation of the output. Here, we derive one of them by assuming ${\lim }_{x\to 0}\frac{\psi \left(x\right)}{x}=c$. We can rewrite (5) in the form,

from which we can easily find an implicit formula for the estimate,

where ω _i,θ is commonly called the weight function. The iterative solution of the M-estimator is usually called W-estimator.

3.1 Robust c-prototypes

As we see in section 2.1, the FCM algorithm is based on a least squares objective function. It is well known that the least square approach is highly sensitive to aberrant points, which implies that FCM gives unsatisfactory results when it is applied to data sets contaminated with noise and outliers. To solve this problem, Frigui and Krisnapuram [16, 17] proposed an algorithm called Robust c-prototypes with the following objective function,

There is one robust loss function ρ _j associated with each cluster. As conventional algorithm, the optimization of this objective function under the same conditions, yields the following updated expressions for the membership matrix and the vector with the centers of the groups:

It is clear that this algorithm can be seen as a generalization of the M-estimator for estimating the prototypes of the FCM algorithm. The expression for the upgrade of the prototypes vector (centers of the groups) is obtained by solving Equation 10, iteratively, which is equivalent to solving the W-estimator (7). By convenience, our proposal makes two significant changes of this algorithm. First, we use the same loss function for all groups. After, the minimization of the objective function with respect to the prototype vector is prevented from entering a new iterative cycle to solve the expression (7) (we must remember that the FCM algorithm is an iterative method), so instead of using the W-estimators, we use the influence function ψ. As was established in section 2.2, the second definition of the M-estimators is based on the existence of the influence function ψ(x, θ) = ∂ ρ(x, θ)/ ∂ θ. Thus, it is not possible to rewrite the expression (10) in terms of the influence function for the W-estimators. Based on these two arguments, the objective function can be written as follows:

In this algorithm and its extension to the color space, we use a gradual assignment of the data to the groups as the conventional FCM algorithm with a value of m = 2. To optimize J _rf with respect to the membership matrix, we use the Lagrange multiplier technique to obtain the following updated equation:

To minimize J _rf with respect to the prototype vector, we fix the membership matrix and set the gradient to zero

3.2 Color robust c-prototypes

The above stated algorithm lets us segment grayscale images, but this idea could be extended on any color space. We extend this proposal to segment color images using different color spaces, but this extension is not immediate. For instance, the RGB color space requires to compute U _f and C _j , 2⁸ ‖ 2⁸ ‖ 2⁸ times, this is equivalent to run the classical FCM algorithm in an image of 4,096 × 4,096 pixels. Since this is computationally intractable in practice, it is necessary to use a color space quantization into q _h bins where h is the channel index [12]. Most of the color spaces use three channels (i.e., RGB, HSV, CIELab, and IJK) to describe each pixel in an image, so that we have to define q ₁, q ₂, and q ₃ and we consider that each color component is divided into the same number of bins q ₁ = q ₂ = q ₃ = q.

In this paper, we propose to apply the proposed color robust c-prototypes (ColorRc P) clustering algorithm in the JK chromatic subspace of the IJK color space. This color space provides that a convenient representation should yield distances and provide independence between chromatic and achromatic components [24].

The best know system for color representation is the RGB (red, green, blue) whose components have values in the interval [0,255]. A color specification system based on luminosity, saturation, and hue form a perceptual representation system that is obtained by some transforms applied to the RGB system [19].

Every color q from the RGB system is described by a vector having three scalar components q = (R, G, B). Let z ₁, z ₂, and z ₃ be three real and positive numbers verifying the equality z ₁ + z ₂ + z ₃ = 1. For q ₁ = (R ₁, G ₁, B ₁) and q ₂ = (R ₂, G ₂, B ₂) the following scalar product is defined in the RGB space as follows:

Let there be the following basis in the RGB space: f ₁ = (1, 1, 1), f ₂ = (1, 0, 0), and f ₃ = (0, 1, −1). Using the Gram-Schmidt procedure, it results the following orthonormal basis: e ₁ = (1, 1, 1), $e_2=\left(\sqrt{\frac{1-z_1}{z_1}},-\sqrt{\frac{z_1}{1-z_1}},-\sqrt{\frac{z_1}{1-z_1}}\right)$, and $e_3=\left(0,\sqrt{\frac{z_3}{z_2\left(1-z_1\right)}},-\sqrt{\frac{z_2}{z_3\left(1-z_1\right)}}\right)$. The coordinates I, J, K for the IJK orthonormal coordinates system in the basis e ₁, e ₂, e ₃ are computed by I = 〈q, e ₁〉 = z ₁ R + z ₂ G + z ₃ B, $=q,e_2=\sqrt{z_1\left(1-z_1\right)}\left(R-\frac{z_{2}G+z_{3}B}{z_2+z_3}\right)$, and $K=〈q,e_3〉=\sqrt{\frac{z_2{z}_3}{1-z_1}}\left(G-B\right),$ where z ₁, z ₂, and z ₃ can be taken as constants with a value of $\frac{1}{3}$ [23]. As a result, the components of the IJK space can be computed as $I=\frac{R+G+B}{3},J=\frac{2R-G-B}{3\sqrt{2}}$, and

In order to avoid the influence of non-uniform illumination, we apply the clustering procedure only for the JK chromatic subspace. By using the color quantification q stated above, we define q ₁ and q ₂ corresponding to the J and K color components. Then, the objective function to minimize is,

subject to ${\displaystyle {\sum }_{q_1=1}^n{\displaystyle {\sum }_{q_2=1}^n{u}_{q_{1}j,q_{2}j}<0}}$, ∀j ∊ {1, …, c} and ${\displaystyle {\sum }_{j=1}^c{u}_{q_{1}j,q_{2}j}<1}$, ∀(q ₁, q ₂) ∊ {1, …, n}².

Taking into account the consideration that each color component J and K are divided into the same number of bins q ₁ = q ₂ = q. We can rewrite the objective function (15) as

Since the gradient of $J_{r_{\mathit{JK}}}$ with respect to u _qj and c _j vanishes when reaching the local optimum, and taking into account the same conditions to minimize the objective function (16), it is easy to show that the optimal updating equations are given by

The proposed ColorRc P clustering algorithm is evaluated here, and its performance has been compared with the following FCM-based segmentation methods: the Hidden Markov random field FCM (HMRFFCM) [9], the histogram thresholding FCM (HTFCM) [10], the fuzzy hyper-prototype clustering (FHCS) [11], the quantized FCM (QFCM_S2) [12], the usual-initialization FCM (UFCM) [13], the color-clustering FCM (CFCM) [14], and the single point iterative weighted FCM (SWFCM) [15]. We also compare our proposal with other image segmentation techniques beyond FCM such as the penalized inverse expectation maximization (PIEM) which extracts the features for each pixel using the Gabor filter, and the classification of pixels in different regions is done by the expectation maximization (EM) algorithm [25] and the segmentation by clustering then labeling (SCLpost) which uses small homogeneous regions to adopt a 3D feature vector obtained through the color space, and then, the clustering is done by a hybrid approach that combines the mean-shift with a semi-supervised discriminant analysis algorithm [26].

When a set of real images are segmented, the criteria used to compare the segmentation performance of various algorithms were the Probabilistic Rand Index (PRI) for evaluation of the results obtained from the tested algorithm to a set of manually segmented images [27], the Variation of Information (VOI) for quantification of the loss of information and the gain between two clusters belonging to the lattice of possible partitions [28], the Global Consistency Error (GCE) for quantification to what extent a segmentation can be viewed as the refinement of the other [29], and the Boundary Displacement Error (BDE) for evaluation of the average displacement error of boundary pixels between two segmented images by computing the distance between the pixel and the closest pixel in the other segmentation [30]. The Kappa coefficient $\widehat{K}$ was used to summarize the segmentation performance when a remote sensing image is segmented [31]. The Kappa statistics considers a measure of the overall accuracy of the image classification and individual category accuracy as a means of actual agreement between classification and observation [32]. Many research works in image segmentation use these indexes to compare the performance of various methods, for this reason, we use these measures to evaluate our proposed algorithms [27–30].

The PRI index is defined as follows [27]:

where N is the number of pixels, S is the segmentation provided by the tested algorithm, c _ij is a Boolean function denoting if l _i ^S = $l_i^{G_k}$, p _ij is the expected value of the Bernoulli distribution for the pixel pair, $l_i^{G_k}$ is the label pixel x _i in the k th manually segmented image, and l _i ^S is the label of pixel x _i in the tested segmentation. The ground-truth set is defined as {G ₁, G ₂, …, G _L }, where L is the number of manually segmented images. The PRI index is in the range [0, 1] where high values indicate a large similarity between the segmented images and the ground-truth; the VOI index [28] is

where $H\left(\Vert \right)=-{\displaystyle {\sum }_{i=1}^c\frac{n_i}{n}\log \frac{n_i}{n}}$ is the entropy, n _i being the number of points belonging to the i th cluster, $I\left(S,G_k\right)={\displaystyle {\sum }_{i=1}^c{\displaystyle {\sum }_{j=1}^c\frac{n_{i,j}}{n}\log \frac{n_i}{n}\frac{n_j}{n}}}$ is the mutual information between two clustering, and n _i,j is the number of points in the intersection of cluster i of S and j of G _k . The VOI is a distance, the smaller the VOI value is, the closer the segmentation obtained and the ground-truth are; the GCE index [29] is

where $E\left(S,G_k,x_i\right)=\frac{\left|R\left(S,x_i\right)\backslash R\left(G_k,x_i\right)\right|}{\left|R\left(S,x_i\right)\right|}$ is a measure of error at each pixel x _i , |^.| is the cardinality, \ is the set difference, and R(S,x _i ) is the set of pixels corresponding to the region in segmentation S that contains the pixel x _i . The better the segmentation S respect to the ground-truth is when the closer GCE is to zero; the BDE index [30] is

where $D_S^{G_k}$ is a distance distribution signature obtained by adding the distances over all points of S; and the Kappa coefficient is computed from [31, 32]

where x _ii are the diagonal entries of the confusion matrix, x _i+ and x _+i indicates the sum of row i and the sum of column i of the confusion matrix, respectively, N _c is the number of elements in the confusion matrix, and ν = 9 is the length of the sample data in a current 3 × 3 window. A value of zero indicates no agreement, while a value of 1.0 shows perfect agreement between the classifier output and the reference data.

Figure 1 depicts a subset of real images from the Berkeley Segmentation Data Set 500 (BSD500) used in the tests. This subset includes a set of real-world color images along with their segmentation maps provided by different people. In these images, a number of difficult aspects come in the image segmentation process, such as natural noise, artifacts, and varying imaging conditions.

Subset of real images from the Berkeley Segmentation Data Set 500 (BSD500). (a) 35070, (b) 42049, (c) 67079, (d) 80099, (e) 113044, (f) 118035, (g) 124084, (h) 135069, (i) 196073, (j) 198023, (k) 208001, and (l) 210088.

Full size image

Table 1 presents the average performance criteria (PRI, VOI, GCE, and BDE) obtained by the proposed ColorRc P algorithm, as well as the comparative algorithms after segmenting all images. We note that the comparative HMRFFCM, HTFCM, and FHCS algorithms are based on the RGB color space [9–11], and the PIEM, QFCM_S2, and SCLpost algorithms are based on the CIELab color space [12, 25, 26]. As can be seen in Table 1, the proposed algorithm has better results than the comparative algorithms in its respective color space.

The optimal parameters for the proposed ColorRc P algorithm are set to initialize randomly, m = 2 for fuzzifier, ϵ = 1E-6 for the tolerance parameter that controls the iteration of the algorithm as well as the quality of the clustering procedure, in our algorithm is the minimum error in the objective function used as a termination criterion, and the number of clusters c are chosen depending on the image to be segmented. The M-estimators used to control the influence of the atypical information are presented in Table 2 [22, 33, 34].

Some visual results are depicted from Figures 2,3,4, and 5. As can be seen in the images, the ColorRc P method can segment the images better and not distort the edges in comparison with other methods.

Segmentations of image 35070 ( c = four regions). (a) HMRFFCM, (b) HTFCM, (c) FHCS, (d) SCLpost (λ = 200), (e) PIEM, (f) QFCM_S2, (g) ColorRc P_Simple-Cut, (h) ColorRc P_Hampel, (i) ColorRc P_{German-McClure}, (j) ColorRc P_Asad, and (k) ColorRc P_Insha.

Full size image

Segmentations of the image 113044 ( c = two regions). (a) HMRFFCM, (b) HTFCM, (c) FHCS, (d) SCLpost (λ = 200), (e) PIEM, (f) QFCM_S2, (g) ColorRc P_Simple-Cut, (h) ColorRc P_Hampel, (i) ColorRc P_{German-McClure}, (j) ColorRc P_Asad, and (k) ColorRc P_Insha.

Full size image

Segmentations of the image 196073 ( c = three regions). (a) HMRFFCM, (b) HTFCM, (c) FHCS, (d) SCLpost (λ = 200), (e) PIEM, (f) QFCM_S2, (g) ColorRc P_Simple-Cut, (h) ColorRc P_Hampel, (i) ColorRc P_{German-McClure}, (j) ColorRc P_Asad, and (k) ColorRc P_Insha.

Full size image

Segmentations of the image 210088 ( c = three regions). (a) HMRFFCM, (b) HTFCM, (c) FHCS, (d) SCLpost (λ = 200), (e) PIEM, (f) QFCM_S2, (g) ColorRc P_Simple-Cut, (h) ColorRc P_Hampel, (i) ColorRc P_{German-McClure}, (j) ColorRc P_Asad, and (k) ColorRc P_Insha.

Full size image

To demonstrate the performance of the proposed ColorRc P clustering scheme in real applications, a sub-region from the remote sensing image ‘Zhalong Nature Reserve’ (see Figure 6a) is segmented [15]. The proposed algorithm is compared with other approaches designed to segment this kind of images, such as the UFCM, CFCM, and SWFCM algorithms [13–15].

Segmentation results on a remote sensing image in c = eight regions. (a) sub-region of Zhalong Nature Reserve on October 21, 2001 [33], (b) UFCM, (c) CFCM, (d) SWFCM, (e) ColorRc P_Simple-Cut, (f) ColorRc P_Hampel, (g) ColorRc P_{German-McClure}, (h) ColorRc P_Asad, and (i) ColorRc P_Insha.

Full size image

As in the first test, the parameters for the ColorRc P algorithm were set as initialized randomly, m = 2 for fuzzifier and a termination criterion of ϵ = 1E-6 based on the minimum error and in the objective function is calculated in each iteration. The remote sensing image was segmented in c = eight regions by using the different M-estimators presented in Table 2.

Table 3 shows the performance segmentation of different algorithms in terms of Kappa coefficient where one can see that the proposed method has a superior value than the other algorithms used as comparative. From Figure 6, the best visual results are provided by the proposed method with different M-estimators with a better preservation of edges and with more homogenous regions.

This paper presents a method to segment color images using the IJK color space; it is based on the robust c-prototypes algorithm, but adapted to use different M-estimators. As could be seen in the quantitative results, the performance of all variants presented offered better results than other algorithms proposed in the recent literature. On the other hand, the visual results depicted that in the segmented images, there are detailed preservation and homogenous regions in comparison with other algorithms.

The authors thank the Instituto Politécnico Nacional de México (National Polytechnic Institute of Mexico) and CONACYT for their financial support.

Authors and Affiliations

1. Instituto Politécnico Nacional, Escuela Superior de Ingeniería Mecánica y Eléctrica, Av. IPN s/n, Edificio Z, acceso 3, 3er piso; SEPI-Electrónica, Col. Lindavista, 07738, Mexico, DF, Mexico

   Dante Mújica-Vargas, Francisco J Gallegos-Funes, Alberto J Rosales-Silva & José de Jesús Rubio

Corresponding author

Correspondence to Francisco J Gallegos-Funes.

Competing interests

The authors declare that they have no competing interests.

Open Access This article is distributed under the terms of the Creative Commons Attribution 2.0 International License ( https://creativecommons.org/licenses/by/2.0 ), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Reprints and permissions