Abstract
Bio-inspired computing has lately demonstrated its usefulness with remarkable contributions to shape detection, optimization, and classification in pattern recognition. Similarly, multithreshold selection has become a critical step for image analysis and computer vision sparking considerable efforts to design an optimal multi-threshold estimator. This paper presents an algorithm for multi-threshold segmentation which is based on the artificial immune systems(AIS) technique, also known as theclonal selection algorithm (CSA). It follows the clonal selection principle (CSP) from the human immune system which basically generates a response according to the relationship between antigens (Ag), that is, patterns to be recognized and antibodies (Ab), that is, possible solutions. In our approach, the 1D histogram of one image is approximated through a Gaussian mixture model whose parameters are calculated through CSA. Each Gaussian function represents a pixel class and therefore a thresholding point. Unlike the expectation-maximization (EM) algorithm, the CSA-based method shows a fast convergence and a low sensitivity to initial conditions. Remarkably, it also improves complex time-consuming computations commonly required by gradient-based methods. Experimental evidence demonstrates a successful automatic multi-threshold selection based on CSA, comparing its performance to the aforementioned well-known algorithms.
1. Introduction
Several image-processing applications aim to detect and classify relevant features which may be later analyzed to perform several high-level tasks. In particular, image segmentation seeks to group pixels within meaningful regions. Commonly, gray levels belonging to the object are substantially different from those featuring the background. Thresholding is thus a simple but effective tool to isolate objects of interest; its applications include several classics such as document image analysis, whose goal is to extract printed characters [1, 2], logos, graphical content, or musical scores; also it is used for map processing which aims to locate lines, legends, and characters [3]. Moreover, it is employed for scene processing, seeking for object detection, marking [4], and for quality inspection of materials [5, 6].
Thresholding selection techniques can be classified into two categories: bilevel and multilevel. In the former, one limit value is chosen to segment an image into two classes: one representing the object and the other one segmenting the background. When distinct objects are depicted within a given scene, multiple threshold values have to be selected for proper segmentation, which is commonly called multilevel thresholding.
A variety of thresholding approaches have been proposed for image segmentation, including conventional methods [7–10] and intelligent techniques [11, 12]. Extending the segmentation algorithms to a multilevel approach may cause some inconveniences: (i) they may have no systematic or analytic solution when the number of classes to be detected increases, and (ii) they may also show a slow convergence and/or high computational cost [11].
In this work, the segmentation algorithm is based on a parametric model holding a probability density function of gray levels which groups a mixture of several Gaussian density functions (Gaussian mixture). Mixtures represent a flexible method of statistical modelling as they are employed in a wide variety of contexts [13]. Gaussian mixture has received considerable attention in the development of segmentation algorithms despite its performance is influenced by the shape of the image histogram and the accuracy of the estimated model parameters [14]. The associated parameters can be calculated considering the expectation-maximization (EM) algorithm [15, 16] or gradient-based methods such as Levenberg-Marquardt, LM [17]. However, EM algorithms are very sensitive to the choice of the initial values [18]; meanwhile, gradient-based methods are computationally expensive and may easily get stuck within local minima [14]. Therefore, some researchers have attempted to develop methods based on modern global optimization algorithms such as the learning automata (LA) [19] and differential evolution algorithm [20]. In this paper, an alternative approach using a bio-inspired optimization algorithm for determining the parameters of a Gaussian mixture is presented.
On the other hand, biological inspired methods can successfully be transferred into novel computational paradigms as shown by the successful development of artificial neural networks, evolutionary algorithms, swarming algorithms, and so on. The human immune system (HIS) is a highly evolved, parallel, and distributed adaptive system [21] that exhibits remarkable abilities that can be imported into important aspects in the field of computation. This emerging field is known as artificial immune system (AIS) [22] which is a computational system fully inspired by the immunology theory and its functions, including principles and models. AISs have recently reached considerable research interest from different communities [23], focusing on several aspects of optimization, pattern recognition, abnormality detection, data analysis, and machine learning. Artificial immune optimization has been successfully applied to tackle numerous challenging optimization problems with remarkable performance in comparison to other classical techniques [24].
Clonal selection algorithm (CSA) [25] is one of the most widely employed AIS approaches. The CSA is a relatively novel evolutionary optimization algorithm which has been built on the basis of the clonal selection principle (CSP) [26] of HIS. The CSP explains the immune response when an antigenic pattern is recognized by a given antibody. In the clonal selection algorithm, the antigen (Ag) represents the problem to be optimized and its constraints, while the antibodies (Abs) are the candidate solutions of the problem. The antibody-antigen affinity indicates as well the matching between the solution and the problem. The algorithm performs the selection of antibodies based on affinity either by matching against an antigen pattern or by evaluating the pattern via an objective function. In mathematical grounds, CSA has the ability of getting out of local minima while simultaneously operating over a pool of points within the search space. It does not use the derivatives or any of its related information as it employs probabilistic transition rules instead of deterministic ones. Despite its simple and straightforward implementation, it has been extensively employed in the literature for solving several kinds of challenging engineering problems [27–29].
In this paper, the segmentation process is considered as an optimization problem approximating the 1D histogram of a given image by means of a Gaussian mixture model. The operation parameters are calculated through the CSA. Each Gaussian contained within the histogram represents a pixel class and therefore belongs to the thresholding points. In order to compare the segmentation results with other optimization methods, the number of elements in the Gaussian mixture (classes) is considered already known or given by the user. The experimental results, presented in this work, demonstrate that CSA exhibits fast convergence, relative low computational cost, and no sensitivity to initial conditions by keeping an acceptable segmentation of the image, that is, a better mixture approximation in comparison to the EM- or gradient-based algorithms.
This paper organizes as follows. Section 2 presents the method following the Gaussian approximation of the histogram. Section 3 provides information about the CSA while Section 4 demonstrates the automatic threshold determination. Section 5 discusses some implementation details. Experimental results for the proposed approach are presented in Section 6, followed by the discussion summarized in Section 7.
2. Gaussian Approximation
Let consider an image holding gray levels whose distribution is displayed within a histogram . In order to simplify the description, the histogram is normalized just as a probability distribution function, yielding where denotes the number of pixels with gray level and being the total number of pixels in the image.
The histogram function can thus be contained into a mix of Gaussian probability functions of the form with being the probability of class , being the probability distribution function of gray-level random variable in class , and being the mean and standard deviation of the ith probability distribution function, and being the number of classes within the image. In addition, the constraint must be satisfied.
The mean square error is used to estimate the parameters , , and ; . For instance, the mean square error between the Gaussian mixture and the experimental histogram function is now defined as follows:
Assuming an n-point histogram as in [13] and being the penalty associated with the constrain . In general, the parameter estimation that minimizes the square error produced by the Gaussian mixture is not a simple problem. A straightforward method is to consider the partial derivatives of the error function to zero by obtaining a set of simultaneous transcendental equations [13]. However, an analytical solution is not always available considering the nonlinear nature of the equation which in turn yields the use of iterative approaches such as gradient-based or maximum likelihood estimation. Unfortunately, such methods may also get easily stuck within local minima.
In the case of other algorithms such as the EM algorithm and the gradient-based methods, the new parameter point lies within a neighbourhood distance of the previous parameter point. However, this is not the case for the CSA which is based on stochastic principles. New operating points are thus determined by a parameter probability function that may yield points lying far away from previous operating points, providing the algorithm with a higher ability to locate and pursue a global minimum.
This paper aims to compare its segmentation results to other optimization methods that have been applied to similar segmentation tasks. Therefore, the number of elements in the Gaussian mixture (classes) is considered as already known or provided by the user. The number of classes, in most cases, represents the number of objects which are contained in the image. However, if the selected number is lower than the object number, it can be assumed that results actually favour the classification of bigger objects yet bearing the expense of ignoring small subjects.
3. Clonal Selection Algorithm
In natural immune systems, only the antibodies (Abs) which are able to recognize the intrusive antigens (nonself cells) are to be selected to proliferate by cloning [21]. Therefore, the fundament of the clonal optimization method is that only capable Abs will proliferate. Particularly, the underlying principles of the CSA are borrowed from the CSP as follows:(i) maintenance of memory cells which are functionally disconnected from repertoire,(ii) selection and cloning of most stimulated Abs,(iii) suppression of nonstimulated cells,(iv) affinity maturation and reselection of clones showing the highest affinities, (v) mutation rate proportional to Abs affinities.
From immunology concepts, an antigen is any substance that forces the immune system to produce antibodies against it. Regarding the CSA systems, the antigen concept refers to the pending optimization problem which focuses on circle detection. In CSA, B cells, T cells, and antigen-specific lymphocytes are generally called antibodies. An antibody is a representation of a candidate solution for an antigen, for example, the prototype circle in this work. A selective mechanism guarantees that those antibodies (solutions) that better recognize the antigen and therefore may elicit the response are to be selected holding long life spans. Therefore, such cells are to be named memory cells .
3.1. Definitions
In order to describe the CSA, the notation includes boldfaced capital letters indicating matrices and boldfaced small letters indicating vectors. Some relevant concepts are also revisited below: (i) antigen: the problem to be optimized and its constraints (circle detection),(ii) antibody: the candidate solutions of the problem (circle candidates),(iii) affinity: the objective function measurement for an antibody (circle matching),
The limited-length character string d is the coding of variable vector x as d = encode(x); and x is called the decoding of antibody d following decode(d).
Set I is called the antibody space; namely, . The antibody population space is thus defined as where the positive integer is the size of antibody population which is an m-dimensional group of antibody d, being a spot within the antibody space I.
3.2. CSA Operators
Based on [30], the CSA implements three different operators: the clonal proliferation operator (), the affinity maturation operator (), and the clonal selection operator (). is the antibody population at time k that represents the set of antibodies a, such as . The evolution process of CSA can be described as follows:
3.2.1. Clonal Proliferation Operator ()
Define where , and is a -dimensional identity column vector. Function round() gets to the least integer bigger than . There are various methods for calculating . In this work, it is calculated as follows: where is called the clonal size. The value of is proportional to the value of . After clonal proliferation, the population becomes where
3.2.2. Affinity Maturation Operator ()
The affinity maturation operation is performed by hypermutation. Random changes are introduced into the antibodies just like it happens in the immune system. Such changes may lead to increase in the affinity. The hypermutation is performed by the operator which is applied to the population as it is obtained by clonal proliferation .
The mutation rate is calculated using the following equation [25]: beingthe mutation rate, being the objective function value of the antibody () as it is normalized between [0,1], and being a fixed step. In [31], it is demonstrated the importance of including the factor into (3.7) to improve the algorithm performance. The way modifies the shape of the mutation rate is shown by Figure 1.
The number of mutations held by a clone with objective function value is equal to , considering as the length of the antibody—22 bits are used in this paper. For the binary encoding, mutation operation can be done as follows: each gene within an antibody may be replaced by its opposite number (i.e., 0-1 or 1-0).
Following the affinity maturation operation, the population becomes where is the operator as it is defined by (3.7) and applied onto the antibody .
3.2.3. Clonal Selection Operator ()
Define, for all , as the antibody with the highest affinity in , then , where is defined as where .
Each step of the CSA may be defined as follows. (1) Initialize randomly a population (Pinit), a set of candidate solutions of subsets of memory cells (M) which is added to the remaining population (), with the total population being , with M holding n memory cells.(2) Select the n best individuals of the population to build A(k), according to the affinity measure (objective function).(3) Reproduce () population A(k) proportionally to their affinity with the antigen and generate a temporary population of clones Y(k). The clone number is an increasing function of the affinity with the antigen (3.1).(4) Mutate () the population Y(k) of clones according to the affinity of the antibody to the antigen (3.4). A maturated antibody population Z(k) is thus generated.(5) Reselect () the best individuals from Z(k) and A(k) to compose a new memory set .(6) Add random novel antibodies (diversity introduction) to the new memory cells M to build .(7) Stop if any criteria are reached; otherwise return to (2.2).
Figure 2 shows the full draw of the CSA. The clone number in Step 3 is defined according to (3.1). Although a unique mutation operator is used in Step 5, the mutated values of individuals are inversely proportional to their fitness by means of (3.7); that is, the more Ab shows a better fitness, the less it may change.
The similarity property [32] within the Abs can also affect the convergence speed of the CSA. The idea of the antibody addition based on the immune network theory is introduced for providing diversity to the newly generated Abs in M, which may be similar to those already in the old memory M. Holding such a diverse Ab pool, the CSA can avoid being trapped into local minima [30], contrasting to well-known genetic algorithms (GAs) which usually tend to bias the whole population of chromosomes towards only the best candidate solution. Therefore, it can effectively handle challenging multimodal optimization tasks [33–36].
The management of population includes a simple and direct searching algorithm for globally optimal multimodal functions. This is also another clear difference in comparison to other evolutionary algorithms, like GA, because it does not require crossover but only cloning and hypermutation of individuals in order to use affinity as selection mechanism. The CSA is adopted in this work to find the parameters , , and of Gaussian functions and their corresponding threshold values for the image.
4. Determination of Thresholding Values
The next step is to determine the optimal threshold values. Considering that the data classes are organized such that , the threshold values are obtained by computing the overall probability error of two adjacent Gaussian functions, yielding considering is the probability of mistakenly classifying the pixels in the (h + 1)th class belonging to the hth class, while is the probability of erroneously classifying the pixels in the hth class belonging to the (h + 1)th class. are the a priori probabilities within the combined probability density function, and is the threshold value between the hth and the (h + 1)th classes. One value is chosen such as the error is minimized. By differentiating with respect to and equating the result to zero, it is possible to use the following equation to define the optimum threshold value : considering
Although the above quadratic equation has two possible solutions, only one of them is feasible; that is, a positive value falling within the interval.
5. Implementation Details
In this work, either an antibody or an antigen will be represented (in binary form) by a bit chain of the form with c representing a point in an L-dimensional space: The clonal selection algorithm can be stated as follows.(1) An original population of N individuals (antibodies) is generated, considering the size of 22 bits.(2) The n best individuals based on the affinity measure are selected. They will represent the memory set.(3) Such n best individuals are cloned m times.(4) Performing a hypermutation of the cloned individuals which follows the proportion inside the affinity between antibodies and antigens and generating one improved antibody population.(5) From the hypermutated population, the individuals with the higher affinity are to be reselected.(6) As for the original population, the individuals with the highest affinity are replaced, improving the overall cells set.
Once the above steps are completed, the process is started again, until one individual showing the highest affinity is found, that is, finding the minimum stated in (2.3). At this work, the algorithm considers 3 Gaussians to represent the same number of classes. However, it can be easily expanded to more classes. Each single Gaussian has the variables (with ) representing the Hamming shape-space by means of an 22 bits word over the following ranges: with g representing the grey level and h is the histogram of the grey level image. Hence, the first step is to generate the initial individual population of the antibody population by means of with representing the bit string assigned to each of the initial individuals . Later, in order to perform the mapping from binary string to real value, we use As to find the corresponding real value for r, with representing .
The population is set to 100 individuals (N), including the best 20 individuals (n). The 20 selected individuals are cloned 10 times (m). The corresponding mutation probability is proportional to the resulting error by (2.3). The algorithm is thus executed until the minimum possible value of (2.3) is reached.
6. Experimental Results
6.1. Presentation of Results
In this section, two experiments are tested to evaluate the performance of the proposed algorithm. The first one considers the well-known image of the “The Cameraman” shown in Figure 3(a) as its corresponding histogram is presented by Figure 3(b). The goal is to segment the image with 3 different pixel classes. During learning, the CSA algorithm adjusts 9 parameters in this test. In this experiment, a population of 100 () individuals is considered. Each candidate holds 9 dimensions, yielding with N representing the individual's number, in this case, 100. The parameters () are randomly initialized.
(a)
(b)
The experiments suggest that, after 130 iterations, the CSA algorithm has converged to the global minimum. Figure 4(a) shows the obtained Gaussian functions (pixel classes), while Figure 4(b) shows the combined graph. The layout in Figure 4(b) suggests an easy combination of the Gaussian functions to approximate the shape of the graph shown in Figure 3(b), representing the histogram of the original image. Figure 5 shows the segmented image whose threshold values are calculated according to (3.4) and (3.5).
(a)
(b)
In order to test the performance, the algorithm gets to analyze the image shown in Figure 6 whose histogram exhibits a remarkable problem (a set of peaks) regarding class approximation. Such image, due to its complexity, is considered as a benchmark image for other algorithms, including classical approaches as in [7, 10] or intelligent algorithms as in [11, 12]. For this particular image, after 190 generations, the algorithm was capable to achieve the minimum approximation value J (see (2.3)), considering three different classes. Figure 7 shows the approximation quality; meanwhile, Figure 8 presents the obtained segmented image.
(a)
(b)
(a)
(b)
The third experiment considers a new image known as “The scene” shown by Figure 9(a). The histogram is presented in Figure 9(b). Now, the algorithm considers 4 pixel classes. The optimization is performed by the CSA algorithm resulting in the Gaussian functions shown by Figure 10(a). Figure 10(b) presents the combined graph from the addition of such Gaussian functions.
(a)
(b)
(a)
(b)
After the optimization by the CSA algorithm, the added layout including all 4 Gaussian functions is obtained as shown by Figure 11(a). It is also evident that the resulting function approaches the original histogram as Figure 11(b) shows the resulting image after applying the segmentation algorithm.
(a)
(b)
6.2. Comparing the CSA versus the EM and LM Methods
In order to enhance the algorithm’s analysis, the proposed approach is compared to the EM algorithm and the Levenberg-Marquardt (LM) method which are commonly employed for computing Gaussian mixture parameters. The comparison focuses on the following issues: sensitivity to initial conditions, singularities, convergence, and computational costs.
6.2.1. Sensitivity to the Initial Conditions
In this experiment, initial values of the mixture model for all methods are randomly set while the same histogram is taken in account for the approximation task. Final parameters representing the Gaussian mixture (considering four different classes) after convergence are reported. All algorithms (EM, LM, and CSA) are executed until no change in parameter values is detected. Figure 12(a) shows the image used in this comparison while Figure 12(b) pictures its histogram. All experiments are conducted several times in order to assure consistency. Table 1 exhibits the parameter values (, , , ) of the obtained Gaussian mixtures, considering the two initial conditions in which the highest contrasts were found. Figure 13 shows the segmented images obtained under such initial conditions. Further analysis on Table 1 shows the acute sensitivity of the EM algorithm to initial conditions. By such sensitivity, it is observed in Figure 13 where a clear pixel misclassification appears in some sections of the image. In case of the LM method, although it does not present a considerable difference in comparison to optimal values, its deviation shows that it is prone to get trapped into a local minimum. On the other hand, the CSA algorithm exhibits the best performance as its parameter values fall the nearest to the optimal approximation performance.
(a)
(b)
6.2.2. Convergence and Computational Cost
The experiment aims to measure the number of iterations and the computing time spent by the EM, the LM, and the CSA required to calculate the parameters of the Gaussian mixture after considering three different benchmark images. Figure 14 shows the images used in the experiment. Such images are selected, since they are employed in the standard segmentation literature. All the experiments consider four classes. Table 2 shows the averaged measurements as they are obtained from 20 experiments. It is evident that the EM is the slowest to converge (iterations), and the LM shows the highest computational cost (time elapsed) because it requires complex Hessian approximations. On the other hand, the CSA shows an acceptable tradeoff between its convergence time and its computational cost. Finally, Figure 14 below shows the segmented images as they are generated by each algorithm. By analyzing the images (a)–(c) in Figure 14, it is clear that the CSA approach presents better results when the segmented images are compared with the original ones. In case of the EM and LM algorithms, several stains are identified in regions where the intensity level must be considered homogenous.
(a)
(b)
(c)
7. Conclusions
In this paper, an automatic image multi-threshold approach based on the clonal selection algorithm (CSA) is proposed. The segmentation process is considered to be similar to an optimization problem. The algorithm approximates the 1-D histogram of a given image using a Gaussian mixture model whose parameters are calculated through the CSA. Each Gaussian function approximating the histogram represents a pixel class and therefore one threshold point.
Experimental evidence shows that the CSA has a compromise between its convergence time and its computational cost when it is compared to the expectation-maximization (EM) method and the Levenberg-Marquardt (LM) algorithm. Additionally, the CSA also exhibits a better performance under certain circumstances (initial conditions) on which it is well reported in the literature [14, 18] that the EM has underperformed. Finally, the results have shown that the stochastic search accomplished by the CSA method shows a consistent performance with no regard of the initial value and still showing a greater chance to reach the global minimum.
Although Table 2 indicates that the CSA method can yield better results in comparison to the EM and gradient-based methods, notice that the aim of our paper is not intended to beat all segmentation algorithms which have been proposed earlier but to show that the artificial immune systems can effectively serve as an attractive alternative to evolutionary algorithms which have been employed before to successfully segment images.