Image Thresholding Improved by Global Optimization Methods

ABSTRACT

Image thresholding is a common segmentation technique with applications in various fields, such as computer vision, pattern recognition, microscopy, remote sensing, and biology. The selection of threshold values for segmenting pixels into foreground and background regions is usually based on subjective assumptions or user judgments under empirical rules or manually determined. This work describes and evaluates six effective threshold selection strategies for image segmentation based on global optimization methods: genetic algorithms, particle swarm, simulated annealing, and pattern search. Experiments are conducted on several images to demonstrate the effectiveness of the proposed methodology.

Introduction

Image segmentation (Bhanu and Lee 2012; Gonzalez and Woods 2010; Rosenfeld 2013; Schwartz and Pedrini 2006; Silva et al. 2008) plays an important role in several image processing and computer vision tasks. Thresholding is one of the simplest techniques for extracting objects in an image from its background. Thresholding techniques (Ali, Ahn, and Pant 2014; Beauchemin 2013; Sezgin 2004; Singla and Patra 2015) can be classified into bilevel and multilevel categories. In bilevel thresholding (Zou et al. 2012; Zou, Liu, and Zhang 2013), a threshold value is determined to segment the image into two regions, one corresponding to the object in the image and another representing the background. In multilevel thresholding (Ali, Ahn, and Pant 2014; Sarkar, Das, and Chaudhuri 2015) more than one threshold value is determined to segment the image into background and objects present in the image.

Several image thresholding approaches (Dirami et al. 2013; Kapur, Sahoo, and Wong 1985; Otsu 1979; Weszka, Nagel, and Rosenfeld 1974; Ye et al. 2012) have been proposed in the literature. Many of them (Otsu 1979; Weszka, Nagel, and Rosenfeld 1974) are based on shape information of the histogram to select threshold values. In ideal cases, the threshold can be selected at the bottom of a valley between two peaks representing objects and background, respectively. However, it is often difficult to detect valleys precisely in real images. Other methods attempt to optimize an objective function, such entropy maximization (Kapur, Sahoo, and Wong 1985), Bayesian error minimization (Ye and Danielsson 1988), Bayesian error minimization (Ye and Danielsson 1988), and intra-class variance minimization (Otsu 1979).

One of the main difficulties related to the use of complex segmentation algorithms is that they have several parameters that are usually defined empirically by users, which affects the quality of segmentation and can vary significantly for different problems.

Since the threshold values defined in a segmentation process can considerably influence the detection and location of objects present in the image, it is important to have a robust method for selecting the best thresholds. In this paper, the parameters used in a number of different image thresholding techniques are automatically defined through four global optimization methods (Floudas 2013; Horst and Tuy 2013).

The main contributions of the work include the definition of an effective fitness function to evaluate the candidate solutions based on an image quality measure and the automatic selection of appropriate image thresholds through optimization approaches.

The remainder of the paper is organized as follows. “Background” section briefly describes concepts and works related to image thresholding techniques. The proposed method for selecting the best parameters in a number of thresholding methods based on optimization methods is presented in “Proposed methodology ” section. Experimental results are described and discussed in “Experimental results” section. “Conclusions” section presents the conclusions and directions for future work.

Background

The main purpose of thresholding (Beauchemin 2013; Dirami et al. 2013; Gonzalez and Woods 2010; Ye et al. 2012) is to extract pixels belonging to objects from an image background. Image thresholding techniques are commonly categorized into global and local approaches (Glasbey 1993; Trier and Jain 1995). In global thresholding, a single threshold value T is used for all the image pixels, which is suitable for cases where the pixel values corresponding to object components and background are consistent along the entire image. In local thresholding, a different threshold T is adaptively selected for each pixel according to local image characteristics, which is more suitable to accommodate non-uniform lighting conditions in the image.

Global thresholding

The simplest thresholding is to select an intensity value as a threshold, such that the pixel values below this threshold become 0 (black) and the values above this threshold become 1 (white). If T is the global threshold of the original image , then

(1)

A widely used global thresholding approach was developed by Otsu (Otsu 1979), which automatically clusters the image pixels into two classes, background and foreground. The method minimizes the intra-class variance (the variance within the class) or, conversely, maximizes the inter-class variance (the variance between the classes).

Considering the range of intensity levels, the intra-class variance can be defined as the weighted sum of the variances of each class

(2)

where and . Value is the variance of the pixels in the background (below the threshold), whereas is the variance of the pixels in the foreground (above the threshold).

The inter-class variance can be defined as

(3)

where and

Local thresholding

Local thresholding techniques are designed to overcome the limitations found in the global thresholding due to changes in lighting conditions or presence of local shadows.

Bernsen’s method (Bernsen 1986) calculates the local threshold value based on the mean value of the minimum and maximum intensities of pixels within a window, defined as

(4)

where and are the minimum and maximum intensity values, respectively, within an region centered at .

Niblack’s method (Niblack 1986) calculates the threshold value for each pixel based on the local mean and standard deviation, expressed as

(5)

where and are local mean and standard deviation, respectively, within an region of pixel . The region size should be small enough to preserve local details of the image, but also large enough to suppress noise. The value k is used to adjust the fraction of total pixels that belong to foreground object.

Sauvola and Pietaksinen’s method (Sauvola and Pietaksinen 2000) calculates the threshold value for each pixel as

(6)

where and are local mean and standard deviation, respectively, within an region of pixel . Sauvola and Pietaksinen suggested the values and .

Phansalskar, More, and Sabale’s method (Phansalskar, More, and Sabale 2011) is variation to Sauvola and Pietaksinen’s method (Sauvola and Pietaksinen 2000) to deal with low-contrast images. The local threshold is calculated as

(7)

where local mean and standard deviation, respectively, within an region of pixel . They suggested the values , , , and .

The mean method (Jain, Kasturi, and Schunck 1995) selects the threshold as the mean of the local intensity distribution in an region. If the intensity of a pixel is above the mean of its neighborhood, then the pixel is set to object; otherwise, it is set as background.

In a similar way, the median method (Jain, Kasturi, and Schunck 1995) selects the threshold as the median of the local intensity distribution in an region. If the intensity of a pixel is above the median of its neighborhood, then the pixel is set to object; otherwise, it is set as background.

Proposed methodology

To find the best parameter values for the image thresholding problem, we model a fitness function to guide the selection of candidate solutions. The Structural Similarity Index (SSIM) (Wang et al. 2004) is a quality measure calculated between two images f and g is expressed as

(8)

where is the mean of f, is the mean of g, is the variance of f, is the variance of g, and is the covariance of f and g. The SSIM measure is in the range of , such that closer it is to 1, the more similar the images f and g are.

The fitness function is defined as

(9)

where is the set of all possible solutions. Therefore, our goal is to maximize the quality measure computed between an input image f and an estimated image g.

In our work, four global optimization approaches (Arora 2015; Chong and Zak 2013) are used to generate the estimated image g by automatically selecting the best parameter values employed in the image thresholding techniques described in the previous section.

Genetic Algorithms (GAs) (Holland 1975; Koza 1992) are heuristic search algorithms based on evolutionary strategies of natural selection, which can be used to generate solutions to search and optimization problems. First, an initial population is created with a number of individuals. Each individual corresponds to a possible solution, that is, it represents a set of possible values for the thresholding parameters. Each gene represents a specific parameter encoded in a bit string. The parameters values for the initial set of solutions are randomly generated. A solution (individual) is evaluated through the fitness function [Equation (8)], where the value resulting from the segmentation is compared to a reference. During the optimization process, the best solutions are selected and new solutions are created from them (reproduction). For each generation, new individuals are created via crossover to compose the next generation. The mutation operator helps prevent local optimum. The evolutionary process continues until a predefined number of generations is reached, such that the segmentation parameters (genes) of the best solution found in the last generation correspond to the final result.

The Particle Swarm Optimization (PSO) (Clerc 2010; Kennedy and Eberhart 1995) algorithm starts with a random population (called swarm) of candidate solutions (called particles), each one having the parameters to be optimized. At each iteration, every particle adjusts its velocity vector and the best known position of the swarm is updated according to the fitness function [Equation (8)]. The algorithm stores and progressively replaces the best parameters of each particle, as well as the particle that best fits the parameters. The process continues until a predefined number of iterations is performed.

The Simulated Annealing (SA) (Brooks and Morgan 1995; Dowsland and Thompson 2012) method models the physical process of heating a material and then slowly decreasing the temperature to reduce defects, thus minimizing the system energy. The algorithm starts with a randomly chosen initial solution. At each iteration, a new point is randomly generated. The distance of the new point from the current point is based on a probability distribution with a scale proportional to the temperature. The algorithm accepts not only new points that produce a better solution, but also points that generate a worse solution. This reheating mechanism helps prevent local optimum, such that other possible solutions can be explored globally. As the algorithm proceeds, a cooling schedule (Kirkpatrick, Gelatt, and Vecchi 1983) is selected to systematically decrease the temperature, reducing the extent of its search to converge to an optimum.

Pattern Search (PS) (Hooke and Jeeves 1961) is a numerical optimization approach that does not require derivatives (gradient) of the fitness function to be optimized. A pattern is a finite set of points (vectors) used for the algorithm to determine which points to search at each iteration. Therefore, the method searches a set of points around the current point to find that one where the value of the fitness function generates a better result than the value at the current point. The algorithm stops when a predefined number of iterations is performed, the variation in the fitness function from one iteration to the next is below a tolerance, or the distance between the point found at one successful poll and the point found at the next one is below a tolerance.

Experimental results

The algorithms proposed in our work were implemented in MATLAB under Linux operating system. Experiments conducted on several images demonstrated the effectiveness of the image thresholding techniques improved by the optimization strategies. Due to space limitations, results for only few images are reported here.

Table 1 reports some results obtained with the optimization methods for selecting the parameters for image thresholding. Values for the computed SSIM measure, percentage of foreground pixels (black), and estimated parameters are shown for each image thresholding algorithm. The size n for the neighborhood region used in the algorithms ranged from 1 to 30.

For the GA, the population size was defined as 50, the crossover probability was 0.8, the mutation probability was 0.05, and the number of generations G was 10. For the PSO, the population size was defined as 50 and the inertia range was set from 0.1 to 1.1. In the SA algorithm, the initial temperature was set to 100. In the PS algorithm, the initial mesh size was set to 1.0 to perform a local search in a neighborhood around the mesh points, whereas the tolerance on variable was defined as 0.000001 to interrupt the iterations if both the change in position and the mesh size are less than such tolerance.

Images generated with the application of the best estimated parameters are shown in Figure 1. It is possible to observe that the fitness function was able to produce adequate qualitative segmentation results. Moreover, the SSIM measures achieved through the image thresholding methods are higher when compared to those obtained with traditional global thresholding Otsu’s method (Otsu 1979).

Figure 1. Images resulting from the application of the best estimated parameters to each thresholding method. The values of SSIM are provided for each image by taking the original one as reference.

It is worth mentioning that the PS optimization method, as it can be observed from Table 1, provided the best qualitative results in relation to the other three methods.

Table 1. Comparison of the results obtained with four different optimization methods.

Image	Method	Result	GA	PSO	SA	PS
Fiducial (640 × 840 pixels)	Bernsen (1986)	SSIM	0.283080	0.283050	0.283050	0.283050
		Foreground (%)	35.91	35.87	35.87	35.87
		N	22	21	22	22
	Niblack (1986)	SSIM	0.288407	0.288657	0.287688	0.288657
		Foreground (%)	38.89	38.98	36.11	38.98
		n	11	10	14	10
		k	–0.000522	–0.0224538	–0.0671209	–0.022460
	Mean (Jain, Kasturi, and Schunck, 1995)	SSIM	0.154173	0.288516	0.288516	0.288516
		Foreground (%)	51.31	39.55	39.55	39.55
		n	7	10	10	10
	Median (Jain, Kasturi, and Schunck, 1995)	SSIM	0.090706	0.092082	0.092082	0.092082
		Foreground (%)	56.24	55.83	55.83	55.83
		n	29	30	30	30
	Sauvola and Pietaksinen (2000)	SSIM	0.289609	0.286833	0.287376	0.290410
		Foreground (%)	36.96	34.91	35.07	32.87
		n	10	18	16	7
		k	0.00406674	0.0142321	0.0125171	0.014707
		R	–33.904	101.883	86.4602	16385
	Phansalskar, More, and Sabale (2011)	SSIM	0.289617	0.289777	0.274469	0.290501
		Foreground (%)	34.27	35.20	23.46	32.98
		n	9	9	7	7
		k	0.0160981	0.0103334	0.178103	–0.000289
		R	49.8369	88.3678	–169.032	1.00513
		p	5.15045	–693.918	16.7733	–0.013382
		q	28.2679	553.543	133.475	0.000457
Watch (768 × 1024 pixels)	Bernsen (1986)	SSIM	0.103089	0.103081	0.103089	0.103089
		Foreground (%)	39.06	39.02	39.06	39.06
		n	17	16	17	17
	Niblack (1986)	SSIM	0.111335	0.111426	0.111416	0.111465
		Foreground (%)	20.41	20.97	20.71	21.00
		n	24	23	22	25
		k	–0.388632	–0.370083	–0.3779	–0.371017
	Mean (Jain, Kasturi, and Schunck, 1995)	SSIM	0.101979	0.106767	0.106767	0.106767
		Foreground (%)	41.54	38.80	38.80	38.80
		n	8	10	10	10
	Median (Jain, Kasturi, and Schunck, 1995)	SSIM	0.0726361	0.0724491	0.0726361	0.0726361
		Foreground (%)	54.85	55.40	54.85	54.85
		n	30	29	30	30
	Sauvola and Pietaksinen (2000)	SSIM	0.125003	0.125006	0.122471	0.125091
		Foreground (%)	20.23	19.90	17.12	19.31
		n	6	7	4	6
		k	0.054655	0.0566141	0.0632998	0.0530624
		R	57.9847	646.569	88.9129	984.746
	Phansalskar, More, and Sabale (2011)	SSIM	0.12487	0.125195	0.118832	0.125426
		Foreground (%)	19.51	19.97	19.84	19.88
		n	6	7	16	6
		k	0.0780287	0.0624377	0.0406067	0.000976562
		R	41.3679	121.909	–10.3443	0.986542
		p	9.16471	9603.67	–53.8454	–0.03578
		q	24.3086	1102.03	189.716	–0.00683975

For each image thresholding, values for the computed SSIM measure, percentage of foreground pixels (black), and evaluated parameters are shown.

Conclusions

The selection of threshold values for image segmentation is typically performed manually based on empirical assumptions; however, this process may significantly affect the quality of the final result.

This work presented and discussed the application of four global optimization methods (GAs, particle swarm, SA, and PS) to six image thresholding techniques.

The automatic threshold values estimated by the algorithms, according to a fitness function based on image quality, were able to generate highly promising results. The proposed method can be used as a preprocessing step in other tasks, such as object location and image interpretation.

As directions for future work, we intend to investigate the improvement of other image processing techniques, such as filtering, registration, and compression.

Funding

The authors are thankful to the São Paulo Research Foundation (FAPESP grant #2015/12228-1) and the Brazilian National Council for Scientific and Technological Development (CNPq grant #305169/2015-7) for their financial support.

Additional information

Funding

The authors are thankful to the São Paulo Research Foundation (FAPESP grant #2015/12228-1) and the Brazilian National Council for Scientific and Technological Development (CNPq grant #305169/2015-7) for their financial support.