An Overview of the Topological Gradient Approach in Image Processing: Advantages and Inconveniences

Image analysis by topological gradient approach is a technique based upon the historic application of the topological asymptotic expansion to crack localization problem from boundary measurements. This paper aims at reviewing this methodology through various applications in image processing; in particular image restoration with edge detection, classification and segmentation problems for both grey level and color images is presented in this work. The numerical experiments show the efficiency of the topological gradient approach for modelling and solving different image analysis problems. However, the topological gradient approach presents a major drawback: the identified edges are not connected and then the results obtained particularly for the segmentation problem can be degraded. To overcome this inconvenience, we propose an alternative solution by combining the topological gradient approach with the watershed technique. The numerical results obtained using the coupled method are very interesting.


Introduction
The goal of topological optimization is to find the optimal decomposition of a given domain in two parts: the optimal design and its complementary. Similarly in image processing, the goal is to split an image in several parts, in particular, in image restoration the detection of edges makes this operation straightforward.
The goal of this work is to show how it is possible to solve some image analysis problems using topological optimization tools. Let us recall that the basic and first idea is based on the work of Amstutz et al. 1 in which the authors propose a solution for crack detection problem using the topological gradient approach. Then, this approach was extended to different problems in image processing like restoration, classification, segmentation, inpainting, and tomography application 2-6 . This paper presents an overview of the applications of the topological gradient approach to major image processing problems with an illustration of both inconveniences and advantages of this technique in image analysis.
Let v be a given noisy image defined in a domain Ω ⊂ Ê 2 and u the restored image. We recall that a classical way to restore the image u is to solve the following partial differential equations PDE problems − div c∇u u v in Ω, ∂ n u 0 on ∂Ω, where c is a small positive constant and n denotes the outward unit normal to ∂Ω. This method is well known to give poor results: it blurs important structures like edges and in order to improve this method, nonlinear isotropic and anisotropic methods were introduced. In our topological optimization approach, c takes only two values: c 0 in the smooth part of the image and a small value ρ on edges. Then classical nonlinear diffusive approaches could be seen as a relaxation of our method that increases the instability of the restoration process. The basic principle of the topological gradient approach is presented briefly in Section 2. The other sections are presented as follows The restoration problem with edge detection is presented in Section 3. Section 4 deals with the classification problem. In Section 5, we present the application of the topological gradient approach to tomography problem. Section 6 is devoted to the segmentation problem; we present in this section the major drawback of the topological gradient approach and propose an alternative to this inconvenience. We end this paper with some concluding remarks. Several numerical tests are also presented for all the applications studied.

Topological Gradient Approach: Basic Principle
Let Ω be an open bounded domain of Ê 2 , and let, j Ω J u Ω be a cost function to be minimized, where u Ω is the solution to a PDE problem defined in Ω. For a small ρ ≥ 0, let Ω ρ Ω \ σ ρ be the perturbed domain by the insertion of a crack σ ρ x 0 ρσ n , where x 0 ∈ Ω, σ n is a straight crack and n a unit vector normal to the crack. The topological sensitivity theory provides an asymptotic expansion of j when ρ tends to zero. It takes the general form where f ρ is an explicit positive function going to zero with ρ and G x 0 is called the topological gradient at point x 0 . It is clear that to minimize the criterion j, we have to insert small holes at points where G is negative. Using this gradient type information, it is possible to build fast algorithms. In most applications, a satisfying approximation of the optimal solution is reached at the first iteration of the optimization process. A topological sensitivity framework allowing to obtain such an expansion for general cost functions has been proposed in the work of Masmoudi 1,7 .  Figure 1: a Initial 512×512 Lena image, b restored image using the linear diffusion method, c restored image using a nonlinear diffusion method, and d restored image using the topological gradient method.
We recall that a classical way to restore an image u from its noisy version v defined in a domain Ê 2 is to solve the following PDE problem where c is a small positive constant, called the conductivity, ∂ n denotes the normal derivative and n is the outward unit normal to ∂Ω. The idea of the topological gradient approach is to consider only two values for the conductivity c: c 0 in the smooth part of the image and 0 or rather undefined in its complementary the edges . For this reason, classical nonlinear diffusive approaches, where c takes all the values of 0, c 0 , could be seen as a relaxation of our method. By enlarging the set of admissible solutions, relaxation increases the instability of the restoration process, and this could explain why our method is so efficient with respect to the algorithm cost.

Application to Image Restoration
The goal of this section is to use the topological gradient approach as a tool for detecting edges for image restoration 5, 6 . In fact, an image can be viewed as a piecewise smooth function and edges can be considered as a set of singularities. The basic idea consists then in applying the topological gradient approach used for crack detection 1 . More precisely, let Ω ρ be the perturbed domain obtained by the insertion of a crack σ ρ . The perturbed solution

3.1
The variational formulation of problem 3.1 is given by where a ρ is the bilinear form, defined on and l ρ is the linear form defined on L 2 Ω ρ by Edge detection is equivalent to look for a subdomain of Ω where the energy is small. So, our goal is to minimize the energy norm outside edges The variation of the cost function j ρ J ρ u ρ induced by the insertion of the small crack σ ρ is given by the topological gradient theory 1 : with Journal of Applied Mathematics 5 where u 0 is the solution of the direct problem 3.1 such that ρ 0 and v 0 is the solution of the following adjoint problems

3.8
We note here that the solutions u 0 and v 0 are calculated in the initial domain without cracks, in the sense that the crack size is equal to zero. Then, the topological gradient can be rewritten in the following ways G x, n M x n, n , 3.9 where M x is the symmetric matrix given by We can deduce that G x, n is minimal when the normal n is the eigenvector associated to the smallest i.e., most negative eigenvalue of the matrix M x . In the following, this eigenvalue will be considered as the topological gradient. Our algorithm consists in inserting small heterogeneities in regions where the topological gradient is smaller than a given threshold α < 0. These regions are the edges ω ρ of the image. The algorithm is as follows.

Restoration Algorithm
i Initialization: c c 0 .
ii Calculation of u 0 and v 0 the solutions of the direct 3.1 and adjoint 3.8 problems respectively.
iii Computation of the 2 × 2 matrix M and its lowest eigenvalue λ min at each point of the domain.
iv Set

3.11
v Calculation of u 1 solution to problem 3.1 with c c 1 .
For numerical experiments, we considered, in a first step, grey level images see Figure 1 . We note that the perturbed image is obtained by adding a gaussian noise σ 20 to the original image, which grey level takes its values in the interval 0, 255 . We present in Figure 1 images restored by different approaches. More precisely, the restored image using the topological gradient approach is presented and compared to both linear and nonlinear diffusion methods 8 . Notice that the convergence of the nonlinear approach is achieved 6 Journal of Applied Mathematics after 53 iterations. Our method can be interpreted as a linear isotropic diffusion method defined in the complementary of edges. In the anisotropic case, not only the amount of diffusion is adapted locally to the data but also the direction of smoothing. Like in the anisotropic nonlinear method, we smooth along image edges while inhibiting smoothing across edges. This can be achieved by replacing the diffusivity scalar by a diffusivity tensor. To show the efficiency of our method, we test an anisotropic algorithm with the topological gradient as follows

Anisotropic Topological Gradient Algorithm
i Initialization: c c 0 .
ii Calculation of u 0 and v 0 : solutions of the direct 3.1 and adjoint 3.8 problems.
iii Computation of the 2 × 2 matrix M and its lowest eigenvalue λ min at each point of the domain.
iv Set c Ac τn A , with A τ, n and c τn c 0 0 0 ε , where τ and n represent the tangent and normal directions.
v Calculation of u 1 solution to problem 3.1 .
In Figure 2, we show that there is no significant difference between the results obtained using the isotropic or the anisotropic approach with the topological gradient method. Our two algorithms are tested on several images. All numerical tests show that the results obtained using the two algorithms are very similar. We show in Figure 2 a zoom of the 512 × 512 Lena image and its restoration by the nonlinear method presented in 8 , and the isotropic and the anisotropic topological gradient algorithms.
We consider now color images such that the RGB model is used. Then, the image I is defined by a 3D vector I x ∈ Ê 3 which represents the intensity of the three colors: red, green, blue, and since each component can be treated as a grey level image, then the topological gradient approach can easily be used, provided the cost function to be minimized has no coupling terms between the several channels. For the restoration problem, one usually set the cost function to the square norm of the gradient of the image. If we consider the most simple case, that is, the standard L 2 norm of the gradient, then the cost function can be easily written as a sum of the three independent terms, one for each channel. Then, the color image problem can be decomposed into three channel-by-channel independent problems. We mention here that the three channels can be supposed to be coupled, and a tensor gradient of the image can be considered as it is proposed in the work of Di Zenzo 9 . Some other coupling techniques have been defined in order to process multichannel images. We can cite for instance the Beltrami diffusion-type process, in which images are considered as manifolds. A nonlinear structure tensor is defined and locally adjusted according to a gradient-type measure and allows one to control the diffusion process 10-12 .
For the numerical tests, we first consider in Figure 3 the natural extension of the topological gradient method to multi-channel images, using the standard gradient norm and we compare this approach with the Ambrosio-Tortorelli version of the Mumford-Shah functional.
Then to improve the restoration process, we propose in Figure 4 to present numerical results using the Di Zenzo approach coupled with the topological gradient approach. image. The noisy image is obtained with an additive gaussian noise, such that the signalto-noise ratio SNR is equal to 10.54. The third image corresponds to the restored image by our algorithm. The SNR of the restored image is 19.23. The identified edges given in the fourth image, correspond to the pixels where the most negative eigenvalue of the matrix M is smaller than a given threshold.
We must note here that, for the case corresponding to the classical topological gradient method, the solution u of 2.2 can be written as u u 1 , u 2 , u 3 , where u k are solutions of the following problems: for k 1, 2, 3, which implies that the cost function j can also be decomposed into three terms, one for each channel, as follows: with 3.14 and the topological gradient G x can be rewritten in the following way: where M x is the symmetric 2 × 2 matrix given by with u k 0 and v k 0 being the solutions of the direct and adjoint problems with respect to each channel k.
On the other hand, for the case corresponding to the topological gradient method using the Di Zenzo gradient, a multispectral tensor associated to the image vector field is where f is the function defined by where M x is the 2 × 2 symmetric matrix given by where H is given by and f is defined by 3.18 .

Application to Image Classification
This section is concerned with the problem of classifying an original image using n predefined classes C i , 1 ≤ i ≤ n, by choosing the color intensity as a classifier. More precisely, let u 0 be the original image defined on an open set Ω of Ê 2 . Our goal is to find n color C i , 1 ≤ i ≤ n, and a partition Ω i , 1 ≤ i ≤ n, of Ω such that u 0 is close to C i in Ω i . The classified image u will then be defined by where {Ω i } i 1,...,n are defined by Ω i x ∈ Ω; x belongs to the ith class .

4.2
We recall that the variational approach consists in defining a cost function measuring the root mean square difference between the original image and the classified image The minimization of J is easy, because for each point x ∈ Ω, we only have to find i x arg min{|u 0 x − C i |; i 1, . . . , n} and add x to subset Ω i x . This can be called the closest class algorithm because each pixel of the original image is assigned in the classified image to its closest class. Inspired by the work of Samson et al. 14 in which the authors propose a classification model coupled with a restoration process, we propose to use the topological gradient approach applied to image restoration problem for both supervised and unsupervised classification problems 3, 4 . Our idea is to still consider the PDE problem given by 2.2 and instead of considering c 0 or c ρ from a numerical point of view on the edge set and c c 0 elsewhere, we consider c ρ on the contours of the image and c 1/ρ elsewhere. As ρ is supposed to be small, then u and u 0 are almost the same.
Otherwise, c 1/ρ, which implies that the PDE defined in 2.2 is nearly equivalent to Å u 0 with Neumann boundary condition normal derivative equal to zero, and this will provide a really smooth image. The classification problem is then as follows.

Supervised Classification Algorithm
i Initialization: c c 0 .
ii Calculation of u 0 and v 0 the solutions of the direct 3.1 and adjoint 3.8 problems respectively.
iii Computation of the 2 × 2 matrix M and its lowest eigenvalue λ min at each point of the domain.
iv Set v Calculation of u 1 solution to problem 3.1 with c c 1 .
vi Application of the k-mean classification algorithm to u 1 .
As in the restoration process, the classification algorithm is applied to color images by decomposing the image in the RGB color space and then we deal separately with the three component images. The topological gradient of the sum can be considered as the sum of the three topological gradients with respect to each channel. Figure 5 shows the classification of a 500 × 380 color image, using 30 color levels whereas the original image has 256 3 different colors , without regularization in a first time, then with regularization. We should mention here that the control of the regularization process can be made by considering progressively different values of c : c c 0 /ρ, c c 0 /ρ 2 , c c 0 /ρ 3 , and so on.
On the other hand, the unsupervised classification corresponds to a classification problem in which the classes are not given. In this case, it is possible to determine them in an optimal way, still by using the topological gradient method. The idea is to study the impact of changing the value of a class C i : C i 1 or C i − 1 on the cost function which implies that for each class, we add or substract, depending on which variation provides the most negative variation one to the value of the class. This algorithm has been applied in the previous numerical tests in order to determine the optimal value of the classes: in Figure 5, the values of 3, 2, and 5 classes with respect to the channels RGB have been determined with this algorithm. One may notice that this algorithm needs at least that the number of classes is known. But in unsupervised classification, it is not the case. The idea is then to add another term in the cost function defined in 4.3 , measuring the number of classes. This can be seen as a regularization term, because we usually do not want to get too many classes. We propose here to rewrite the cost function given by 4.3 as follows: where α is a positive regularization coefficient. It is then possible to bound n as the two terms of the cost function are positive. The unsupervised classification algorithm is as follows

Unsupervised Classification Algorithm
i Initialization with n 1.
ii Determination of the optimal value.
iii While n is smaller than its upper bound: iv set n : n 1 and determine the optimal values of the new n classes.
v Stop if the value of the cost function at the optimum with respect to C i is larger than the previous optimal value with n − 1 classes.
This algorithm has been applied in Figure 5 in order to find the optimal number of classes and their optimal values.

Application to Tomography Problem
In this section, we use the topological gradient as a tool for the reconstruction problem in tomography. We recall that a standard approach for regularizing the ill-posed problem of tomographic imaging consists in the following optimization problem:  Figure 6: Evolution of the PSNR in dB of the reconstructed image by the topological gradient method, as a function of the level of the regularization coefficient c 0 .
where A denotes a system matrix defining the discrete Radon transform, g is the measured projection, φ is a regularization functional, and λ represents a parameter which controls the tradeoff between a good fit to the data and the smoothness of the solution. Inspired by the work of. Chan et al. 15 , we propose to consider the following minimization problem: where c is a positive constant and α is equal to 1 on the edges and to 2 elsewhere. Using the same idea of the topological gradient approach, we consider, for a given regular function g, the following direct problem:

5.3
where n denotes the outward unit normal to ∂Ω ρ and A and A # are respectively the discrete Radon transform and the dual discrete Radon transform.  Let v be the solution of the adjoint problem ∂ n v 0 on ∂Ω.

5.4
By considering the same algorithm applied for image restoration, we compute the reconstruction solution by solving the following proble: ∂ n f 0 on ∂Ω, where c 0 is a given constant. In order to obtain the optimal value of the coefficient c 0 , we have studied the behavior of the SNR of the reconstructed image as a function of c 0 . In fact, we have experimentally confirmed that the quality of the reconstructed image is first increasing rapidly with c 0 and then decreasing slightly, as we can see it in Figure 6. The optimal value is c 0 30, and it has then been considered for our numerical tests. We have tested our method on the well known Shepp-Logan head phantom, for which projections data have been computed using a discrete Radon transform. Figure 7 a shows the original image. A gaussian noise is added to the data with a signal to noise ratio SNR 22, the sinogram is represented in Figure 7 b . Figure 7 c shows the reconstruction result PSNR 26.18 . We also present in this paper a numerical comparison using several methods. Table 1 shows the PSNR peak signal to noise ratio, in dB , SNR signal-to-noise ratio, in dB , SSIM structural SIMilarity , and MSE mean squared error of the reconstructed images for several methods: FBP filtered backprojection , FBP Hamming with a Hamming filter on the FBP image , and the topological gradient approach using either the TV total variation norm or both TV and L 2 norms. The sinogram has a noise level such that the SNR is 24.5. We can see in Table 1 that the best results for all indicators are obtained with the hybrid L 1 /L 2 norm within the topological gradient approach. The peak SNR is higher of more than 3 decibels, the structural similarity index is larger than 94% and the mean squared error is smaller than 1%. Note that both topological gradient methods provide much better results than the filtered backprojection, which is not surprising.

Application to Image Segmentation
We propose in this part a new algorithm for the segmentation problem based on the topological gradient approach. Through the different and various applications which were presented, we can easily conclude that the topological gradient approach is a very interesting method in image analysis. First, the topological gradient method can easily be applied to many problems in image analysis. Second, we should mention that we obtain excellent numerical results and the computing time is very interesting and confirms the theoretical complexity O n · log n of our algorithms, where n is the size of the image. However, the topological gradient approach presents a major drawback: the identified edges are not connected and then the results obtained for the segmentation problem can be degraded, particularly for complex images. So, the main idea in this part is to take advantage of the topological gradient efficiency, to detect the main contours with an interesting computational cost the topological gradient algorithms require only three system resolutions and to find an alternative solution for the discontinuous edges detected by the topological gradient approach. The approach proposed in this paper is based on mathematical morphology 16, 17 and more precisely on the watershed transformation 18 . A major advantage of the method proposed is that it attenuates the oversegmentation usually observed with the watershed technique. This attenuation is due to the fact that the topological gradient provides a global analysis of the image then the almost unwanted contours due to the noise added to a given image are effectively reduced by our approach. The basic idea observed using this coupled method is the use of the topological gradient instead of the morphological gradient classically used with watersheds. The property of the topological gradient to be less sensitive than the morphological or classical gradients to noise and small variations of the image, will contribute significantly to give the best identification of the edges and then to attenuate the oversegmentation usually obtained using watershed techniques. We illustrate in this work the efficiency of the combined method proposed for the segmentation problem. Other numerical tests and comparisons with other methods can be found in 19 . In this part, we consider two real images which have various difficulties curves, circles, straight lines, etc. . The first image is a 302 × 280 fruit-basket gray level image and the second one is a 399 × 246 gray level image which represents a road scene. Figure 8 shows the segmentation result of the two previous images, according to the coupled method for image segmentation based on a watershed algorithm 18 and the topological gradient approach previously presented in this paper. The results shown in We can clearly observe the attenuation of the obtained regions such that the main edges of the two images are conserved. One should remark that we can either define the watershed of the function u or of its gradient: the difference between the two definitions is that in the first case we obtain the influence zones of the processed image, while the second case gives the image edges. It is well-known that in both cases, the watershed gives an oversegmentation, and to avoid this drawback, the alternative to combine it with the topological gradient approach seems to be a good choice to overcome this inconvenience as shown in Figure 8.
Finally, due to a large number of approaches used for the segmentation problem, it clearly appears that it is important to compare our experimental results with methods already proposed in the literature. Particularly, we propose to compare our method with an active contour model based on the level set approach, which is well known to be an efficient method, extensively used in many applications over the last decade. The basic idea in active contour models is to evolve a curve subject to constraints from a given image, in order to detect different objects in that image. To achieve this goal, we start with a curve around the object to be detected, the curve moves toward its interior normal and has to stop on the boundary of the object. This was the first idea of classical snakes and active contour models proposed by 20 . Our numerical tests are based on a level set method proposed by 21 and given by the following evolution equation: ∂u ∂t g |∇I| div ∇u |∇u| α |∇u| ∇g, ∇u , in 0, ∞ × Ω, 6.1 with boundary and initial conditions given by ∂u/∂n 0 on 0, ∞ × ∂Ω and Φ 0, x, y Φ 0 x, y in Ω, α ≥ 0 is a positive given coefficient, g |∇I| is an edge detector function defined in our numerical tests by g s 1/ 1 s 2 , and φ 0 represents the initial level set function. We refer the reader to 8 for more details about geodesic active contours and level set methods. Figure 9 shows the numerical results obtained. We should mention that in order to obtain the final contour illustrated in Figure 9 e , the algorithm requires 500 iterations and a computational time of 78 seconds. We present in Figure 9 some solutions as time evolves: during the evolution, the initial curve is shrinking and stopping as soon as it is close to an object boundary and spilling in order to detect the others objects.

Conclusion
We presented in this paper many applications in image processing of the topological gradient approach with an illustration of both advantages and inconveniences of this technique. We proved that the topological gradient approach provides an excellent frame for solving different image processing problems, same if an alternative way to overcome its major drawback to give not connected contours is to combine it with other approaches. It has been applied to image restoration, edge detection, image classification, and image segmentation for both grey level and color images. It has also been applied to tomography problem, and our main goal consists in applying it to other real life problems mainly that in all the applications proposed, we obtained good results and the computing time is very short. We propose in a forthcoming paper to perform the application of the topological gradient approach to image segmentation, eventually by combining it with other approaches and particularly using marker criteria. In fact, as the edges are detected in regions where the topological gradient is the most negative, then it suffices to extract some points which belong to the edge set and then use these points as a selected marker set. We also intend to extend this work to color image segmentation and three-dimensional segmentation.