A novel PDE-based image denoising approach is proposed in this paper. One designs here a nonlinear filter for image noise reduction based on the diffusion flow generated by the porous media equation , where is a nonlinear continuous function of the form , . With respect to standard 2D Gaussian smoothing and some nonlinear PDE-based filters, this one is more efficient to remove noise from degraded images and also to reduce “staircasing” effects and preserve the image edges.

1. Introduction

Image noise removal constitutes a very important process, often representing the first step of image analysis [1, 2]. Image denoising with feature preservation represents still a focus in the image processing domain, remaining a serious challenge for researchers. A robust denoising technique has to successfully reduce the amount of noise while preserving the edges and has no blurring effect on the processed image.

The most common noise results from the image acquisition system can be modeled as Gaussian random noise in most cases. Gaussian noise represents statistical noise having the probability density function equal to that of the normal distribution. Numerous image denoising models have been introduced in the past few decades. The conventional filters such as averaging filter, median filter, and 2D Gaussian filter are efficient in smoothing the noise, but also have the disadvantage of blurring image edges [1, 2]. For this reason, some nonlinear image filters, which produce more satisfactory noise reduction results and preserve better the integrity of edges and detail information, were introduced [3].

The nonlinear noise removal techniques based on Partial Differential Equations (PDEs) have been extensively studied in the last two decades [46]. The PDE models are the best candidates for a very efficient image denoising. We proposed several robust PDE-based models for noise removal and image restoration in our previous works [79]. In this paper we develop a novel PDE-based image denoising approach based on nonlinear diffusion.

The idea of using the diffusion in image processing arose from the use of the Gaussian filter in multiscale image analysis. Nonlinear diffusion methods reduce noise and enhance contours in images [811]. Numerous nonlinear PDE denoising approaches based on diffusion have been introduced since the early work of Perona and Malik in 1987 [10]. The anisotropic diffusion developed by them, also called Perona-Malik diffusion, was intended to smooth the image while preserving its edges. Many denoising schemes based on this influential work have been proposed in the last 25 years [10, 11]. More closely of our approach is perhaps Kacur and Mikula scheme [12] which we will briefly present later on.

In the next section of this paper we describe the porous media diffusion filtering model. Then, the explicit version of the fast diffusion filtering scheme is presented in Section 3. Our denoising technique is compared with Perona-Malik scheme and Kacur-Mikula scheme in the fourth section. Other method comparisons and the numerical experiments performed using the proposed approach are discussed in Section 5. The paper ends with a section of conclusions.

2. The Porous Media Diffusion Filtering

Numerous PDE image denoising and segmentation approaches are based on the following nonlinear diffusion equation [10, 11]: where is the original noised image and is the filtered image. Here is the image domain, is the image boundary, and is a given continuous and monotone function. Very often, is a gradient mapping of , where is a convex continuous function.

To (1) boundary conditions of Dirichlet type: on or Neumann: on are taken. The special case, with ,  , leads to the so-called total variation model. There is a large variety of nonlinear PDE models inspired by (1) which attempt to eliminate the staircasing effect or to recover images with edges and the total variation model is one of the most efficient. Here we propose a nonlinear filter based on the flow generated by the porous media equation:

The function can be multivalued but maximal monotone, that is, for all and the range of is all of . This equation is well known as mathematical model for fluid diffusion in porous media and the typical example is where and . The case models low diffusion processes and fast diffusion. We note that the model (1) was already used in our previous works for restoring (inpainting) the grayscale images [7, 8].

Given the domain , we denote by , , , , the standard spaces of Lebesgue integrable functions and the Sobolev spaces on . We set ,  , and the dual space of with the norm denoting . Consider here (2) with Dirichlet homogeneous boundary conditions, that is, with , and on , under the following assumption: is monotonically increasing and .

Here represents a bounded domain of with a sufficiently smooth boundary, . By strong solution to (3) on , we mean an absolutely continuous function such that the strong derivative in exists a.e. on ,   a.e. and where , and the Laplace operator is considered in the sense of distributions on domain . It should be said that a strong solution to (4) is also a weak solution in the classical sense. Let us consider the nonlinear operator defined by with and rewrite (4) as the infinite dimensional Cauchy problem:

We note [13, 14] that the nonlinear operator is maximal monotone (equivalently, m-accretive) in and that it is just the subdifferential of the function

in the sense of convex analysis. Here, . As a matter of fact, the distribution space is the basic functional space to study problem (4). Then, we have the following.

Theorem 1. Let . Then, (4) (equivalently, (6)) has a unique strong solution satisfying Moreover, is given by the exponential formula

It should be recalled (see [13]) that this existence result extends to nonhomogeneous Cauchy problem ,  .

Theorem 1 is a consequence of the standard existence theory for the infinite dimensional Cauchy problem associated with maximal monotone (m-accretive) operators in Hilbert spaces [13, 14]. It should be emphasized that Theorem 1 applies to (1) as well as to other nonlinear diffusion equations of monotone type arising in the filtering theory. In particular, (9) amounts to saying that the finite difference scheme, where , is convergent to . In the linear case, that is, , this scheme is at origin of high-order Sobolev gradient method developed in [15, 16]. Formula (9) is known in the literature as the Crandall-Liggett exponential formula.

The flow is a continuous semigroup of contractions in in this case or in , in other situation, and the denoising effect of (6) is due to smoothing effect of on initial data, a property which is a characteristic to semigroups generated by the nonlinear operators of the gradient type , where is a convex and lower semicontinuous function. Since Theorem 1 applies to any initial data , which is a quite general space of distributions on , the filter can be used for smoothing very noised images .

It should be emphasized that, for filtering purpose, it suffices to compute at any time since the smoothing effect does not increase with . On the contrary, since for we have [13] where is a steady-state solution to A, that is, , one might expect that, for a large t, the filtered image could be far away from the original image . This means that in the algorithm (10), the best filtered image can be obtained for ,  , and arbitrarily small. If one wants to restore the image via step descent algorithm, then (6) should be replaced by or, in its discrete version, because, by (11), , where , that is,

This procedure was used in [7] for restoring and inpainting significantly degraded images . We also have the generation formula . This means that the operator- is the infinitesimal generator of the semigroup [12]. The smoothing performance of the filter and its standard deviation is made precisely below.

Theorem 2. Let be such that . Then If , then

Proof. Formula (15) is immediate by the continuity of the function from to . As regards formula (16), we have ,  , and this yields ,  , while by equation ,  , we have, by multiplying with in , , for all .

We also note that for (4) the maximum principle is also valid. More precisely, we have the following.

Theorem 3. Assume that , a.e. in . Then,   .  .

Proof. We multiply (4) by and integrate on to get , as claimed.

For denoising procedure, the filter , generated by (2) in the fast diffusion case, that is, for ,   (like the classical Gaussian filter), has the disadvantage that it is not localized as in this case the solution to (2) propagates with infinite speed. Moreover, has extinction in a finite time which depends on the norm of . The situation is different in the low diffusion case ,  , where the propagation is with finite speed, and so the denoising flow is localized.

In (2) one might consider the Neumann boundary condition on . In this case the operator is given by (5) where and Theorem 1 remains valid in the present situation; the basic space is instead of .

3. The Explicit Scheme of the PDE Algorithm

Now, let us consider the explicit version of the scheme (10), obtained by the fast diffusion model. Namely, ,  . We take , where and , and . Thus, the finite difference scheme leads to the following iterative process: where ,  , and .

Using the iterative scheme (17) with some properly selected parameters , , and , respectively, the initial noised image is successfully filtered in steps. Therefore, the obtained represents the final image enhancement result. The choice of in our simulation was dictated by the numerous tests we have performed in specific examples.

A proper selection of the parameter values is very important and cannot be a priori defined. It turns out that the selection of a very large number of iterations, for example, using a K > 40 value, could produce a blurring effect on the processed image, while considering a very small value, such as K < 5, could result in an unsatisfactory image denoising result. A great value increases also the computational complexity of this filtering process, producing a much higher computation time.

Also, using a large enough value, such as , could increase the degradation of the noised image. A very small parameter, such as , produces no visible smoothing results. Obviously, the parameter must satisfy the condition , for a successfully noise removal.

4. Comparison with Other Noise Removal Techniques

The denoising model discussed above is generated by the minimization problem given by (14), where and is the potential function corresponding to nonlinear diffusivity function . The function is a penalty term which forces the restored image to stay close to the initial image .

We compared our nonlinear diffusion-based noise removal technique with Perona-Malik denoising scheme [5, 10, 17] and other more general techniques (see [12]). These approaches based on (1) reduces to the minimization problem:

The fact that the distance from to is taken in the norm which is considerably weaker than the -norm considered in the Perona-Malik model as well as in most of the denoising models [911] has the advantage that it allows to work with very degraded initial images which practically are not represented by Lebesgue integrable functions but by distributions. However, it should be said that our model has a considerable better smoothing effect than that proposed by Perona and Malik. Indeed, as seen above, the solution to minimization problem satisfies the nonlinear elliptic equation: with boundary conditions (or Neumann flux boundary conditions as the basic space is , the dual of , instead of . As seen earlier, the minimization problem can be replaced by the evolution equations (4) or (13)).

This means that which indicates a smoothing effect comparable with the standard one, but for a suitable choice of one can avoid the “staircasing” effect which is common in denoising procedure with high smoothing effect.

A reformulation of problem (14) and, implicitly, of the dynamic model (3) allows a closer comparison of this model with that inspired by the Perona-Malik classical denoising technique [10]. If we denote by the operator ,   and recall that is an isometry from to for , we can rewrite the minimization problem as where is the square root of . The Euler-Lagrange optimality condition for variational problem (20) can be expressed as Keeping in mind that is defined by we may interpret (21) (resp., (20)) in an appropriate sense (1). From this perspective, the denoising model (19) is close to the “total variation model” (1) if is taken as Indeed, in this case, (21) reduces to

However, the present model is more convenient that the “total variation model,” which is constructed in a nonenergetic space (the space of function with bounded variation) and so hard to treat from the computational point of view. As a matter of fact, by the regularization necessary to construct a viable numerical scheme, the “total variation model” loses most of the theoretical advantages regarding the sharp edge detection and elimination of the staircasing effect.

The simulation developed in the next section confirms this important advantage (see the results depicted in Figures 1 and 2). We must emphasize also that there are numerous empirical denoising schemes obtained by modifying the original Perona-Malik model [5, 10, 11, 1719], most of them not well posed from the mathematical point of view. The Kačur and Mikula denoising model [12] is based on the boundary value problem: with Neumann boundary condition on .

Here and are continuous and monotonically increasing, and are Lipschitzian, and is a smoothing kernel. In the special case ,  , and b, independent of problem (25) reduces to (2) with Neumann boundary conditions and replaced by . In general, (25) can be rewritten as with appropriate Neumann boundary conditions. It should be said that in this general form the Cauchy problem for (25) or (26) is not well posed and so the finite difference scheme (10) is not convergent. The model considered here is mathematically well posed. It should be mentioned however that in [12] a semi-implicit scheme is designed which is convergent to a weak solution to (25).

5. Numerical Experiments and Method Comparisons

The PDE-based image denoising technique proposed here has been tested on various image datasets, satisfactory filtering results being obtained. We have filtered hundreds of grayscale images affected by various levels of Gaussian noise using the described approach. An important advantage of our filtering technique is that it performs a robust noise reduction while preserving the image edges [20, 21].

We have obtained the best denoising results for the following set of parameters of the diffusion model provided by formula (18): , which corresponds to the physical model of diffusion in plasma, and . Some grayscale image reduction examples based on these parameter values are provided in Figures 1 and 2, respectively.

Numerous method comparisons have also been performed. The denoising performance of our method has been compared with performances of other noise removal techniques, such as the 2D Gaussian filter, the averaging filter, and Perona-Malik anisotropic diffusion scheme. From the performed numerical experiments it is obvious that the approach introduced here works better than these well-known filtering algorithms, as one can observe in the figures.

In order to assess the performance levels of each image denoising technique, one computes the norm of the error image. Thus, for any initial image and its restored version u, having the size, the error value is calculated. One can see the values of this norm of the image error parameter in the next tables.

In Figure 1(a), the standard image of Lena is displayed in the grayscale form. Then, it is corrupted by an amount of Gaussian noise characterized by parameters 0.2 (mean) and 0.02 (variance), the noising result being displayed in Figure 1(b).

In Figure 1(c) the image denoising result produced by the classic Gaussian 2D filter kernel is displayed, while the noise reduction obtained with an averaging filter kernel is represented in Figure 1(d). The noise removal produced by the Perona-Malik approach is displayed in Figure 1(e), while the denoising result provided by the nonlinear PDE model proposed here is represented in Figure 1(f).

The values of the norm of the error image corresponding to all these image filtering methods are displayed in Table 1. One can see that our anisotropic diffusion-based noise reduction approach performs better than the other filters, minimizing the respective error.

Another image denoising method comparison example is displayed in Figure 2. The same noise removal approaches (standard Gaussian, average filtering, Perona-Malik, and our technique) are applied on the standard Baboon image corrupted by the same amount of Gaussian noise, their results being represented in Figures 2(c)2(f). The corresponding values of the norm of the error image parameter are displayed in Table 2. One can observe the minimum error value obtained by our proposed AD-based noise removal algorithm.

Our filtering algorithm increases the image quality not only by performing a robust Gaussian noise reduction, but also by enhancing the grayscale image contrast. Also, the proposed denoising technique runs quite fast, a digital image being processed in less than 1 s. The time complexity decreases proportionally with the size of parameter.

6. Conclusions

A novel PDE-based image noise reduction technique has been proposed in this paper. We have designed a robust nonlinear image filter based on the flow generated by the porous media equation.

The provided nonlinear diffusion-based denoising technique performs successfully for grayscale images corrupted by 2D Gaussian noise, producing also an improved image contrast. Our method is also an edge preserving noise removal technique, which represents an important advantage.

It has been compared with some influential anisotropic diffusion approaches, like that proposed by Perona and Malik [5, 10], the obtained comparison results being very encouraging for us. Also, it provides better smoothing results than many other denoising techniques, while executing almost as fast as them, given its quite low-time complexity.

We have also tested this PDE-based model on other types of image noise and obtained mixed filtering results. While our technique performs somewhat well for images affected by some types of noise such as speckle [9], it does not succeed in removing properly some other types of noise, like salt and pepper. Therefore, designing some robust filters for other types of image noise will represent the focus of our future work in the image processing domain.

Because of its edge-preserving and contrast-enhancing character, the technique described in this paper can be successfully applied in image analysis domains requiring image object emphasizing and detection.


This work was mainly supported by the project PN-II-ID-PCE-2011-3-0027-160/5.10.2011 financed by UEFSCDI, Romania. It was supported also by the Institute of Computer Science, The Romanian Academy, Iasi, Romania.