Abstract

The Perona-Malik (PM) model is used successfully in image processing to eliminate noise while preserving edges; however, this model has a major drawback: it tends to make the image look blocky. This work proposes to modify the PM model by introducing the Caputo-Fabrizio fractional gradient inside the diffusivity function. Experiments with natural images show that our model can suppress efficiently the blocky effect. Also, our model has good performance in visual quality, high peak signal-to-noise ratio (PSNR), and lower value of mean absolute error (MAE) and mean square error (MSE).

1. Introduction and Some Basic Definitions

Image processing based on partial differential equations (PDEs) is mainly used for smoothing and restoration purposes. Typical PDE techniques for image smoothing regard the original image as initial states of a parabolic process and extract filtered versions from its temporal evolution. In this paper, represents the image intensity values in the position for a time , is the smooth boundary, and is the original image. It is a classical result that for any bounded , the linear diffusion process possesses the solution

For more details, see [1]. However, can be a smooth kernel or a blur kernel. In 1960, Gabor remarked that the difference between the original and the blurred image is roughly proportional to its Laplacian [2]. To formalize this remark, we have to notice that is spatially concentrated and that we may introduce a scale parameter for , namely [2], then

Here, . When gets smaller, the blurring process looks more and more like the heat equation [2]. The linear diffusion equation (1) does not only smooth noise, but it also blurs important features such as edges and, thus, makes them harder to identify. To remove the noise while preserving the edges at best, Perona and Malik proposed the following nonlinear diffusion method (called by them anisotropic diffusion) [3]: where is the original image intensity function, is the smoothed image intensity function in time , denotes the gradient, is the divergence operator, denotes the derivative in normal direction to the boundary, and is the diffusion coefficient function, also called edge-stopping function. The edge-stopping function is a nonnegative decreasing function satisfying two conditions. One is that as , so that the rate of diffusion is high within uniform or inner regions, and the other one is as , so that the diffusion is totally zero across boundaries. The important property of edge functions is that they should have an insignificant value for those gradients that correspond to edges [3]. In the literature, different edge-stopping functions have been proposed. For example, in [3], Perona-Malik used where and is the gradient magnitude threshold parameter that decides the amount of diffusion that take place (see for instance a way to estimate it in [4]). Other expressions for the diffusion coefficient can be found in [1, 5, 6].

The flow function defined as represents the sum of the brightness flow that is generated. The maximum flow is generated at locations where [7].

The Perona-Malik equation (6) and the large amount of its modifications [5, 811] have demonstrated to be able to achieve a good trade-off between noise removal and edge preservation. Unfortunately, edge-stopping functions lead to backward-forward problems that are ill-posed [1]; for this reason, in [8, 12], authors introduced a modification in the diffusion coefficient to obtain a regularized version as follows: where is a bounded domain of with an appropriately smooth boundary, denotes the unit outer normal to , and . The term is the regularized version of , and it is the gradient of the solution at time of the linear heat equation (1). Despite the ill-posedness, the practical results seem to be quite good in most cases. As reported by Weickert in [1], the reason is that the numerical scheme used by Perona and Malik does not correspond to their equation but rather to a time-regularized one which is well-posed this time. Numerically, the mainly observable instability is the so-called staircase (blocky) effect, where a sigmoid edge evolves into piecewise linear segments that are separated by jumps. This effect is visually unpleasant and is likely to cause a computer vision system to falsely recognize as edges the boundaries of different blocks that belong to the same smooth area in the original image [13].

To avoid staircase effects while achieving a good trade-off between noise removal and edge preservation, many authors have proposed to use fourth-order partial differential equations [13]. For example, You and Kaveh proposed in [1317] a fourth-order PDE for noise removal and an algorithm that can remove speckle effect efficiently. In [15], the authors studied the fourth-order telegraph-diffusion equation for image restoration in which the numerical method not only preserves the edge but also avoids the staircase effects, and the existence, uniqueness, and stability were stabilized. The experiments show that all these fourth-order models can improve the peak signal-to-noise ratio, preserve texture, and eliminate the staircase (blocky) effect efficiently (see [18] and the references therein). The problem with the use of fourth-order equations is that they tend to leave the image with isolated black and white speckles (so-called speckle effect). This effect is characterized as pixels whose intensity values are either much larger or much smaller than those of the neighboring pixels [13, 19]. To eliminate the undesirable speckle effect, in [20], the energy functional associated with the Perona-Malik equation is redefined as an increasing function of the absolute value of the image intensity fractional derivative function (defined in the Fourier domain). The corresponding Euler-Lagrange equation is then a fractional-order anisotropic diffusion equation. The proposed pseudo-PDEs will lead to an interpolation between Perona-Malik equations and fourth-order anisotropic diffusion equations in [13]. In [18], a fractional-order Perona-Malik diffusion (FOPMD) equation is proposed. In this model, the integer-order derivative concerning spatial variables of the Perona-Malik diffusion is replaced with the fractional-order Grunwald-Letnikov derivatives (Polubny). The FOPMD model is given by where the fractional order is () and is the smooth gray scale image at time . The fractional-order gradient vector with order is defined as where represents the partial-order derivative of with respect to the variable whose order is . Although model (15) has reported good performance of preserving edges and suppressing staircase and speckle effect, the resulting images still have some artifacts. Therefore, in [21], a new diffusion model named regularized fully spatial fractional-order Perona-Malik diffusion (RFS-FOPMD) is proposed, given by

As observed, equation (19) is obtained by substituting in equation (12), the integer spatial derivative for fractional-order derivative according to Grunwald-Letnikov. The fractional-order is taken as . With this method, results obtained in [18] are improved, partly by improving the performance of edges locating by regularization [21]. Despite this, there are still errors in locating the position of the edges. In [22], it is said that the basic reason for locating edge positions falsely is that the fractional-order gradient module cannot be used as an edge indicator; therefore, the authors in [22] adopted a so-called external fractional-order gradient vector Perona-Malik diffusion by only replacing integer-order derivatives of the external gradient vector to their fractional-order counterparts while keeping first-order derivatives for diffusion coefficient. The model is the same as [9] except for the derivative used in the diffusion coefficient. A novel fractional PDE model is given by Guidotti and Longo in [23]; they address the well-posedness of the following fourth-order model for noise removal, using fractional derivatives defined by Fourier transform. with periodic and . Model (22) is a modification of the You and Kaveh model [13]. As numerical approximation of the equation, they used a scheme based on the Krylov subspace spectral method.

There are many definitions of fractional derivatives (three popular definitions were given by Grunwald-Letnikov (G-L), Riemann-Liouville (R-L), and Caputo). These have been used in numerous fields of science such as study of the anomalous diffusion phenomenon [2426], circuit theory [2729], and image processing [30, 31], among other applications [11, 3248]. Given the discussion above, we consider that using anisotropic diffusion models to eliminate noise in an image, preserving both strong and weak edges and without phenomena such as staircase, speckle, or any type of artifact, is a subject where much remains to be investigated. In image processing, integer-order differentiation operators are often used in edge detection, especially the first order for the gradient (e.g., Roberts, Prewitt, and Sobel) and second order for the Laplacian (e.g., Laplacian of Gaussian). However, the first-order derivative methods generally cause thicker edges, resulting in the loss of image details. The second-order derivative methods have a stronger response to fine detail, but they are more sensitive to noise [49]. To solve this problem, the fractional-order derivative has been introduced in the edge detection methods, with the capability to preserve more low-frequency contour features in the smooth areas, maintain high-frequency marginal features, and also enhance medium-frequency texture details [49, 50]. Many fractional-order operators are used for edge detection, such as the fractional-order Sobel operator [31], fractional-order CRONE operator [51], and quaternion fractional differential operator [52].

Inspired by these ideas, and taking into account that the use of fractional derivatives inside the diffusivity function has proven to be robust in the presence of noise as said in [53], we present a modification to the Perona-Malik model. Instead of the diffusion coefficient controlling the diffusion from the threshold of the integer-order gradient, it does so according to the fractional-order gradient threshold. Our goal by incorporating the fractional gradient into the diffusion coefficient is to avoid undesirable artifacts such as blocky effect while removing noise and preserving edges. To the best of our knowledge, from the consulted literature, this model using the recently defined Caputo-Fabrizio derivative has not been used. Caputo and Fabrizio introduced a new fractional derivative in [54], intending to describe structures with different scales.

We propose to insert the Caputo-Fabrizio fractional gradient into the diffusivity function as follows:

The flow function corresponding to equation (25) is a function depending on two variables and takes the following form: where is given by (9) and and . This means that the graph of function (28) is a surface. For simplicity reasons, the integer-order gradient will be fixed by a real positive value while the fractional-order gradient will be taken as the independent variable. The simulation of function (28) with different values of produces pictures related to Figure 1. From these simulations, we observe that for , the flow function takes higher values. Otherwise, when takes great values, we see a different behavior.

1.1. Some Basic Definitions in Fractional Calculus

Let , such that and . We denote by the space of all continuous functions on with compact support. We denote as

The Sobolev space is defined by

We set

Three popular definitions of fractional calculus were given by Grunwald-Letnikov (G-L), Riemann-Liouville (R-L), and Caputo. Of these, G-L and R-L are the most popular definitions used in digital image processing. Recently, Caputo and Fabrizio introduced a new fractional derivative in [54]. Considering the spatial variable, it is defined as with and . For more details, see [4, 5459].

The remainder of this paper is organized as follows: Section 2 presents discretization of numerical schemes. Experimental results are provided in Section 3. Finally, conclusions are summarized in Section 4.

2. Numerical Schemes

The main purpose of this section is to give an explicit finite difference scheme of the proposed model. For it, we will start discretizing the Caputo-Fabrizio fractional derivative. Next, we will recall the numerical approximation for the PM equation.

2.1. Discretization of the Caputo-Fabrizio Derivative

To solve numerically the new fractional Perona-Malik model, we propose to discretize the Caputo-Fabrizio derivative based on the forward finite difference scheme in the interval (analogously ). Let us take a partition of nodes of the interval , with step and ; then, we obtain

As in digital gray image , the shortest distance of 2-D image on and coordinates is one pixel; then, we put , and from (32), we obtain where

An analogue expression is obtained for .

2.2. Perona-Malik Scheme

Here, we recall the numerical approximation for the PM equation. Let and be the space steps such that , and denoting the time step as , the discretization of the PM equation is given in [13], as where with for the numerical scheme to be stable.

2.3. Numerical Method of the Proposed Model

As we said above, we propose to introduce the Caputo-Fabrizio fractional gradient inside the diffusivity function of the Perona-Malik model where is used in this work. Considering that then the discretization of , using the forward finite difference scheme, is given by

The discrete fractional-order gradient is an 8-dimensional vector, as in [18]: where represents the transpose of the vector and each component is defined as where is given by (34). The discretization of the usual gradient, following the idea in [3], in eight directions, is given by where , , , and are defined as in (37), (38), (39), and (40), respectively, and

The subscripts N, S, E, W, NE, SW, NW, and SE denote the eight directions North, South, East, West, North-East, South-West, North-West, and South-East. Hence, the explicit discretization scheme for equation (42) is given by where

3. Experimental Results

In this section, experimental results are obtained applying the proposed fractional-order model on the images which were corrupted by the Gaussian noise with mean zero and variances and . Gaussian noise is common in images acquired by cameras and telescopes, and it modifies all the pixels of the image. Furthermore, this noise type is one of the most studied in the literature on account of the frequency with which it occurs in practical situations. Our model is compared with the classical Perona-Malik (PM) anisotropic diffusion model. We set and used and . To implement equation (53) numerically, in (47), we considered , which means that in the approximation of fractional derivative, we take into account two points in each direction for experimental calculus. In the literature, there are some methodologies for the estimation of the contrast parameter [61]. However, for reasons of simplification, we will use fixed contrast parameters which are . The number of iterations used is 40. The performance of the two PDE filters has been assessed by using some well-known quality measures such as the PSNR (peak signal-to-noise ratio), MSE (mean square error), and MAE (mean absolute error) which are defined by [60] where and denote the width and height of the image, respectively, is the pixel in the filtered image, and is the pixel in the original image. denotes norm and denotes norm (Euclidean distance). In the comparison of statistic parameters, it is important to note that the larger the PSNR value, the better the statistical result. On the contrary, the smaller the MSE and MAE, the better the result of the model.

Figures 27 show a comparison of denoised images obtained by our proposal and the Perona-Malik models. Tables 13 show the performance of our proposal compared to the Perona-Malik model using the quality measures MAE, MSE, and PSNR. All these tables have been obtained by using . From Tables 13, it can be seen that the PSNR values of the proposed algorithm are higher than those of the Perona-Malik model; thus, it can efficiently delete the noise while preserving the original content of the image. Meanwhile, the values of MAE and MSE of the proposed algorithm are lower, which shows that the proposed method can better preserve texture details in the denoised image and has a better visual effects. Images of the third and fourth rows of Figures 27 have been obtained by using and , respectively. After employing the proposed filter with the two previous mentioned values combined with the above indicated threshold parameters , we observed that the parameter with the best statistic results is when , as shown in Tables 13.

4. Conclusion

In this paper, we proposed a new version of Perona-Malik diffusion (PMD) using the Caputo-Fabrizio fractional gradient inside the diffusivity function. Experiments showed that filtered images by the proposed method look better than results with the Perona-Malik filter. Our technique has demonstrated good performance in visual quality, higher peak signal-to-noise ratio (PSNR), and lower value of mean absolute error (MAE) and mean square error (MSE) than the Perona-Malik filter. Images of the third and fourth rows from Figures 27 have been obtained by using and , respectively. Experiments show that the best statistics are obtained for the fractional derivative of order when , as can be seen in Tables 13. As future work, we want to introduce the convolution between the Caputo-Fabrizio fractional gradient and the Gaussian filter inside the edge-stopping function.

Data Availability

To support this study, the USC-SIPI Image Database was used.

Conflicts of Interest

The authors declare that they have no conflict of interest.

Acknowledgments

This work is supported by Universidad Nacional de Guinea Ecuatorial (UNGE) and Havana University.