Abstract and Applied Analysis

Volume 2014 (2014), Article ID 948392, 16 pages

http://dx.doi.org/10.1155/2014/948392

Research Article

A Total Variation Model Based on the Strictly Convex Modification for Image Denoising

Boying Wu,¹ Elisha Achieng Ogada,1,2 Jiebao Sun,¹ and Zhichang Guo¹

¹Department of Mathematics, Harbin Institute of Technology, Harbin 150001, China

²Department of Mathematics, Egerton University, P.O. Box 536-20115, Egerton, Kenya

Received 18 April 2014; Revised 16 May 2014; Accepted 26 May 2014; Published 18 June 2014

Academic Editor: Adrian Petrusel

Copyright © 2014 Boying Wu et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

We propose a strictly convex functional in which the regular term consists of the total variation term and an adaptive logarithm based convex modification term. We prove the existence and uniqueness of the minimizer for the proposed variational problem. The existence, uniqueness, and long-time behavior of the solution of the associated evolution system is also established. Finally, we present experimental results to illustrate the effectiveness of the model in noise reduction, and a comparison is made in relation to the more classical methods of the traditional total variation (TV), the Perona-Malik (PM), and the more recent D-α-PM method. Additional distinction from the other methods is that the parameters, for manual manipulation, in the proposed algorithm are reduced to basically only one.

1. Introduction

Noise removal, edge detection, contrast enhancement, inpainting, and segmentation have been the subject of intense mathematical image analysis and processing research for nearly three decades. Several methods have been pursued over the passage of time. These include wavelet transform [1, 2], curvelet shrinkage methods [3–6], and variational partial differential equation (PDE) based methods [7, 8]. These methods generate processes that can easily be divided into either linear and nonlinear processes or isotropic and nonisotropic processes [9].

Due its ability to preserve crucial image features, such as edges, nonlinear anisotropic diffusion is favored over isotropic diffusion [10]. Much interest, therefore, has focused on understanding operations and mathematical properties of the nonlinear anisotropic diffusion and associated variational formulations [11, 12], formulation of well-posed and stable equations [8, 11], extending and modifying anisotropic diffusion [11, 13, 14], and studying the relationships that exist between the various image processing techniques [11, 15, 16].

The objective of any image denoising process should not focus only on the removal of noise, but it should also ensure that no spurious details are created on the restored image and that the edges are preserved or sharpened [7, 17, 18]. It is, therefore, necessary to develop formulations which are sensitive to the local image structure, especially edges/contours [19]. Consequently, a number of edge indicators have been proposed and logically grafted into the partial differential equation (PDE) based evolution equations [7, 11, 20].

Some of these PDEs originate from variational problems. For instance, Rudin et al. [8] proposed a minimization functional, widely referred to as the total variation (TV) functional, of the form

where is the fidelity parameter, denotes the noise image, and is an open bounded subset of . TV functionals are defined in the space of functions of bounded variation (BV) and, therefore, do not necessarily require image functions to be continuous and smooth. This fact makes them allow jumps or discontinuities, and hence they are able to preserve edges.

The original TV formulation has certain weaknesses. Firstly, the formulation is susceptible to backward diffusion since it is not strictly convex. Secondly, the numerical implementation cannot be accomplished without additional small perturbation quantity, say , at the denominator [21–23]. Otherwise a spike/singularity is suddenly generated when in the homogeneous regions. This perturbation phenomenon is believed to contribute to some loss of accuracy in the results of the restoration process. Moreover, the additional parameter unnecessarily increases the number of parameters, thereby making it difficult to determine which permutation of parameter values will give optimal result.

Additionally, given that the method is more efficient in preserving edges of uniform and small curvature, it may excessively smoothen and possibly destroy small scale features having more pronounced curvature edges [24]. TV regularization approach may also result in a loss of contrast and geometry of the final output images, even in noise free observed images [9, 25]. Furthermore, TV regularization has difficulties recovering texture, and there is also evidence of enhanced noise when the fidelity parameter is chosen so that texture is not removed [26]. Lastly, the formulation favors piecewise constant solutions. This has the material effect of causing staircases (false edges) on the resultant image, especially from images severely degraded by noise [25].

However, given the strength of TV based techniques, especially in edge preservation, various modifications have been proposed. For instance, Vogel in [22, 27] proposed a total variation penalty method of the form

where , is a linear operator, and is the penalty parameter. The formulation becomes total variation formulation of the Rudin et al. form [8] when and therefore largely suffers the same shortcoming of the ordinary TV formulation. Hence, it has no significant practical advantage over TV.

Strong and Chan [28] proposed an adaptive total variation based regularization model of the form

where is a control factor which controls the speed of diffusion depending on whether the region is homogeneous or an edge. This model demonstrated fairly good results. However, since it is also not strictly convex, it is still susceptible to backward diffusion, which has the potential of introducing blurs in the restored image.

Chambolle and Lions [29] proposed to minimize a combination of total variation and the integral of the squared norm of the gradient and thus have

where is a parameter. In this formulation signals edges, while signals homogeneous regions. This model is successful in restoring images where homogeneous regions are separated by distinct edges but may become sensitive to the thresholding parameter in the event of nonuniform image intensities or heavy degradation [24].

A variable exponent adaptive model was proposed by Chen et al. in [24]. It has the form

where , with proposed as , is the Gaussian filter, and , are fixed parameters. It is observed that the method exploits the benefits of Gaussian smoothing when and the strength of TV regularization when . This method also demonstrates good results and is indeed a great improvement of the earlier models. However, the convolution of the image with a Gaussian kernel before the evolution process diminishes the accuracy of these models. This is because of the introduction of the scale variance of the Gaussian; the scale variance is itself an additional parameter that is subject to manual manipulation [12]. Furthermore, the diffusion process becomes ill-posed if scale variance is too small, while image features become smeared if the Gaussian variance is too large [30]. Therefore, optimal selection of the scale variance remains a challenge. Parameter permutation giving optimal result is a challenge due to the number of such parameters involved.

Further literature surveys attest to the fact that research in effective regularization functionals which have the ability to generate diffusion processes that restore images, while simultaneously preserving critical images features, the analysis, and practical implementation of such models, is still an extremely active concern.

Consequently, in this paper, we propose a new adaptive total variation (TV) formulation for image denoising, which is strictly convex. The only parameter , which is a thresholding parameter to be tweaked, is such that it depends on the evolution parameter and is therefore not as grossly limited as the choice of numerical values. The other parameter is dynamically obtained and is therefore not manually tuned. With all this the number of parameters is reduced to basically only one. Our method approximates TV model, for higher values of the thresholding parameter . Experimental results, presented here, demonstrate the effectiveness of our model over the classical models of TV, PM, and D--PM.

The structure of this paper is as follows. In Section 2, we present the proposed model (6), its properties, justifications of the model, and the associated evolution equation. In Section 3, we give certain preliminary definitions we rely on from time to time in this paper; we also prove the existence and uniqueness of the solution to the minimization problem (6). In Section 4, we discuss the associated evolution equation to the minimization problem (6). Also, we define the weak solution to the evolution problem (12)–(14), derive the estimates for the solution of an approximating problem (32), prove the existence and uniqueness of the weak solution of the evolution problem (12)–(14), and discuss the asymptotic behaviour of the weak solution as . In Section 5, we give the numerical schemes and experimental results to demonstrate the strength and effectiveness of our method. We have also presented a brief comparative discussion of our results with the PM, original TV methods, and D--PM methods. A brief summary concludes the paper in Section 6.

2. Proposed Model

In this section, motivated by the methods of [11, 21, 24], we propose a modified version of total variation denoising model. The model is based on minimization of a total variation functional with a logarithm-based strictly convex modification. First, we present the model and certain important properties of it.

2.1. The New Energy Functional

The new strictly convex energy functional is given as follows:

where

where is a positive parameter, and is the noise image.

Calculating directly, we can obtain the following proposition.

Proposition 1. The function satisfies the following properties:

(i), , and therefore is a strictly nondecreasing function for ;(ii), and therefore is strictly convex;(iii) as ;(iv), , , and therefore is the linear growth.

Compared with the work in [31], the new model satisfies the linear growth condition, which is much more natural. And the existence result for the flow associated with the minimization problem in [31] is only in one and two dimensions because the methods employed there use general results on maximal monotone operators and evolution operators in Hilbert spaces.

Remark 2. In this method, satisfies some minimal hypotheses; namely:

(H1) is a strictly convex, nondecreasing function from to , with ;(H2)because in the homogeneous regions weak gradients should be smoothed, there should be weak penalization in these region. Thus as , ; this signals isotropic diffusion;(H3)because edges, signaled by regions of strong gradients, should be protected, regularization near the edges will be penalized strongly. Thus, as , , ; this not only demonstrates the linear growth nature of the model but also signals the TV edge preservation behaviour of the model.

In [27], the authors obtained

where . Actually, the perturbation is always used to the eliminate singularity of the term in the numerical experiments. With the inclusion or addition of , what is then implemented is not actually the TV formulation, but an approximation of it. And, so, because of the approximation, the time step in the discrete scheme becomes limited only to smaller values, as diminishes [32]. The proposed model beats all these challenges and therefore has an easier design for numerical implementation without the need for lifting parameters. Moreover, it is noticed that, for ,

which is the form of the TV model. This means that any merits accruing from the TV model could be still obtained as becomes larger and larger.

Remark 3. (1) Since TV potential (i.e., ) is only convex (i.e., ), it gives a local minimum, but it cannot guarantee uniqueness of the minimum energy. The proposed energy potential (i.e., ), however, is strictly convex (i.e., ). And, with additional strict convexity in the fidelity part, , the resultant energy functional is strictly convex in both and , thereby giving a global minimum energy and hence guaranteeing uniqueness of the results.

(2) The proposed method seeks to reduce the number of parameters subject to manual manipulation. In the implementation, is made to depend on time, thereby leaving time (evolution parameter) as the only parameter to be tweaked.

Other modifications of the TV model are only convex in and therefore do not guarantee uniqueness of solutions; hence uniqueness of results is not easy to assure. Besides, the evolution system of the models contains component, which runs into a spike when , in smooth regions. This then makes only the implementation of the approximations of the corresponding formulations possible (see [13, 20, 24, 27, 29]).

2.2. The Associated Evolution Equation

The Euler-Lagrange equation for the energy functional (6) is

with

To compute a solution of (10) numerically, it is embedded into a dynamical scheme, where is used as the evolution parameter. Hence, corresponding to the proposed minimization model (6), we have the following evolution system:

Note 1. The modified formulation also does not have the disadvantage of running into a singularity, since the denominator does not become zero. This particular observation makes it more convenient to design a numerical scheme for the model.

Furthermore, to see more clearly the action of the diffusion operator (kernel), we decompose the divergence term on the premise of the local image structures. That is, we break it into the tangential () and normal () directions to the isophote lines. Consequently we have

where

The divergence term is visibly a weighted sum of the directional derivatives, along the normal and tangential directions. Since we are also seeking to preserve the edges and other important features of the image, smoothing in the tangential direction should be encouraged more than that in the normal direction, in the regions near the boundaries (edges). Observe in the above formulation that, as increases, the coefficient of diminishes faster, reducing diffusion in the normal direction, thus preserving edges. The edges are signalled by higher values of the magnitude of gradient of . However, as decreases, there is relatively uniform diffusion in both and directions, thereby achieving isotropic diffusion in homogeneous regions. The homogeneous regions are signaled by diminishing values of the magnitude of gradient of image . The proposed model is, therefore, sensitive to the local image structure. However, TV based denoising model and most other modifications do not have that component. This leads to a situation where the homogeneous parts are processed into piecewise constant regions, whose boundaries reflect staircases or false edges in the image. This situation, for instance, affects the results of a modification by Chen et al. [24] when , and the model becomes the traditional TV model.

3. The Minimization Problem

In this section, we prove the existence and uniqueness of the minimization problem (6). But, first, we give the following preliminary presentations which will guide our reasoning in the subsequent sections of this paper.

3.1. Preliminaries

Definition 4. Let be an open subset of . A function has a bounded variation in if

where denotes the space of such functions. Then -norm is given by

Definition 5 (see [33]). If , then

and is a radon measure, where is the density of the absolutely continuous part of with respect to the Lebesgue measure, , and is the singular part.

Lemma 6 (see [34]). Let be convex, even, nondecreasing on with linear growth at infinity. Also let be the recession function of defined by

Then for and setting we have

This implies that is lower semicontinuous for the weak* topology.

3.2. Existence and Uniqueness of Solution to the Minimization Problem

The linear growth condition makes it natural to consider a solution in the space

However, this space is not reflexive. But sequences bounded in are also bounded in and are therefore compact for weak* topology. Moreover, due to the fact that -space is separable, weak* topology allows us to obtain compactness results even if the space is not reflexive [34–37]. Hence, by denoting weak* topology simply as , we therefore seek the solution to the minimization problem (6) in the space .

Theorem 7. The minimization problem (refer to (6)) has a unique solution .

Proof . From the linear growth property of we have

which implies that

Thus is coercive. Therefore, let be a minimizing sequence in . Then

where denotes a generic constant that may differ from line to line. It clearly follows from above that

This implies that is bounded in . Hence there exists a subsequence of and a function such that

And, by Lemma 6 and the weak lower semicontinuity of the -norm, we have that

Hence, there exists a solution of the minimization problem. Uniqueness of the solution follows from the strict convexity of the functional.

4. Existence and Uniqueness for the Evolution Equation

In this section, we define a weak solution for the evolution equations (12)–(14). That is, if (6) has a minimum point , then it should formally satisfy the Euler-Lagrange equation (12), subject to the boundary conditions (13)-(14). We define a weak solution that we are seeking for the equation system (12)–(14). We then propose an approximating problem and derive certain a priori estimates to the approximating problem. These will aid the process of proving the existence and uniqueness of the weak solution. Finally, we show the large time behavior (asymptotic stability) of the weak solution.

4.1. Definition of Weak Solution

Let , , , and , where and is a solution to (12)–(14). Multiplying (12) by and integrating over , we have

Then applying the convexity condition of , that is, , on both the second term on the left side and on the right side of the above equation gives

Integrating (30) with respect to we get

Now, since , we may choose , where . We observe that the left-hand side of (31) has a minimum at . This shows that is a solution of (12) in the distributional sense. The facts raised above motivate the following definition for a weak solution to (12)–(14).

Definition 8. A function is called a weak solution of (12)–(14) if , , on and satisfies (31) for every and .

4.2. Approximating Energy Functional

Let be an open bounded subset of and let . Then let us consider the following approximating energy functional for :

where

with the following properties:

(1) and ;(2). Hence is strictly convex in .

Remark 9. Since we can have a sequence with

And as a consequence of Theorem  2.9 in [20] we have

From property (iii) of it is easy to see that

This implies that

Hence

On the strength of (32) and Remark 9 let us consider the approximate evolution problem

Let the following bound also hold for the solution of (39)–(41). Then, the following lemma indicates the boundedness of the solution of the above approximate evolution problem.

Lemma 10. Suppose and that is a solution to the problem (39)–(41); then

Proof. By a method similar to that of Zhou in [33], let and define such that

then multiplying (39) and integrating over yields

Observe that the last two integrals in the above equation are nonnegative, based on the definition of and the fact that . Therefore, we have

which indicates that is a decreasing function in . Since we have that

Hence . Conversely, multiplying (39) by , where , and employing a similar argument, we obtain that , . Therefore, , .

Another estimate (bound) is obtained via the following lemma.

Lemma 11. The approximate evolution problem (39)–(41) has a unique solution with such that

for any and with . Moreover

Proof. Since is a lower semicontinuous and strictly convex function, then, by definitions provided in [34, 38, 39], , which is a subdifferential of , is maximal monotone. The approximate problem (39)–(41) has a unique solution such that . Now, to show that is the weak solution to the approximating problem (39)–(41), as governed by (47), we multiply (39) by for with and integrate over and . Thus

Applying convexity condition on both sides of the above equation we obtain

which implies (47).

Then, multiplying (39) by and integrating over and we obtain

which yields

From the above equation together with (41) we obtain

which by (40) and (41) produces (48).

4.3. Existence and Uniqueness of the Solution of the Evolution Problem

Next, using the a priori estimates obtained above, through the approximate problem, we proceed to prove the existence and uniqueness of the solution of the evolution problem (12)–(14).

Theorem 12. Let . Then there exists a unique weak solution , , and such that

where , .

Proof. From Lemmas 10 and 11 we see that

where is some constant that may vary from line to line. Let by Lemma 10 be a bounded sequence of solutions of the problem (39)–(41). Then there exists a subsequence denoted by of and a function , with , such that, as ,

Indeed, from (55), there is a sequence and a function with such that (57) and (58) hold. Observe also that for any , as ,

which indicates that for each

By Lemmas 10 and 11, for each , is a bounded sequence in . Combining that fact with (62) we obtain that for each as

Moreover, (60) follows from the fact that

From (64) is equicontinuous, and from (55) and (58) we have that

It then follows by standard argument that we can have in uniformly in , giving us (59). From Lemma 10 and (58) we obtain that with .

Next we show that, for all , and and for each

To show this end, from Lemma 11, we obtain

From (56) and (58) we deduce that there exists a subsequence such that

Using (57), (59), and (68) and we let ; in (67) we obtain

But we define from the lower semicontinuity theorem that

Using (70) in (69) and integrating over we obtain

which confirms (66).

Now, to complete the proof of the existence of solution to (12)–(14), it remains to pass to the limit as in (71). Replacing by in (48), letting , and using (57)–(59) and (70) we obtain

From Lemma 10 and also from the above inequality we see that is uniformly bounded in and is uniformly bounded in . This means we can extract a subsequence of such that as we have

Now, replacing with in (71), letting , and applying the lower semicontinuity in (70) we obtain

for all , and and for each . Hence is a weak solution to (12)–(14). Replacing by in (72), letting , and using (57)–(59) and (70) we obtain (54).

Uniqueness of the Weak Solution. With reference to the definition of solution inequality as given in (31), let and be two weak solutions to the problem (12) to (14) such that . Then we have, for two solutions, two inequalities:

Now, adding the two inequalities (75) we obtain a more compact inequality given by

This implies that

which implies for a.e. . This confirms the uniqueness of the weak solution.

4.4. Large Time Behavior

Finally, we will assess the asymptotic limit of the weak solution as . The essence of this part is to demonstrate that over time the solution of the evolution problem (12)–(14) ultimately converges to the unique minimizer of the functional (6).

Theorem 13. As the weak solution of the evolution equation (12)–(14) converges weakly in to a minimizer of the functional in (6).

Proof. Let into (31) to obtain

The equation above simplifies to

By taking , and have for each , . Given that is uniformly bounded in , we conclude that there exists sequence and its subsequence is still denoted by such that in and weak* in , as , respectively. Hence dividing inequality (79) by and taking the limit as we obtain

This demonstrates that indeed is a weak solution to (12)–(14) and also the unique minimizer of problem (6).

5. Numerical Experiments

In this section, we present the performance of our method in denoising images involving a Gaussian white noise. We have then compared our results with the ones obtained by the classical methods of PM method [7], TV method [8], and the more recent D--PM method [11].

5.1. Numerical Scheme

In the following two subsections, two numerical discrete schemes, the PM scheme (PMS) and the additive operator splitting (AOS) scheme, have been proposed.

5.1.1. PM Scheme

Here, we have proposed a numerical scheme similar to the original PM method, whereby (12)–(14) are discretized as follows:

and is dynamically determined according to the following discretization scheme:

Hence from (81) and (82) we have

where

for and .

5.1.2. AOS Scheme

In this part, using a AOS scheme, the problem (12)–(14) has been discretized as follows:

where ,

where

where is the set of the two neighbors of pixel (boundary pixels have only one neighbor).

It is observed that AOS schemes with large time steps still reveal average grey value invariance, stability based on extremum principle, Lyapunov functionals, and convergence to a constant steady state [10]. The AOS scheme is less than twice the typical effort needed for the PM scheme, a rather low price for gaining absolute stability [10]. It is worth noting that if we use to approximate (TV kernel), in numerical scheme, the AOS scheme may be unstable because of the small number . Our approximation, however, effectively avoids instabilities arising from such a scenario, as there is no need for lifting the denominator, since it cannot not be equal to zero. This has been an additional motivation for considering AOS for our numerical experiments. The possibility of occurrence of a zero denominator in the evolution problem is a phenomenon that makes the numerical implementation of TV problematic.

5.2. Comparison with Other Methods

The experiments in this work were performed on a Compaq610 computer, having Intel Core 2 Duo CPU T5870 each 2.00 GHz, physical RAM of 4.00 GB, and Professional Windows 8 64-bit Operating System, on MATLAB R2013b. The image restoration performance was measured in terms of the peak-signal-to-noise ratio (PSNR), mean absolute deviation/error (MAE), structural similarity index measure (SSIM), the measure of similarity of edges (PSNRE), and visual effects. The iteration stopping mechanism was based on the maximal PSNR or the lowest MAE. The PSNR and MAE values were obtained using the formula by Durand et al. [40] and are given by

Also, taking the edge map from the evolution equation, we have

And corresponding edge similarity measure is given, according to the formulation by Guo et al. [11], by

The SSIM measures have been obtained according to the formula by Wang et al. [41] given by

where denotes the noise-free image, is the denoised image, is the dimension of image, and yields the gray scale range of the original image. In addition, compares the two images’ mean luminances and . The maximal value of , if and , measures the closeness of contrast of the two images and . Contrast is determined in terms of standard deviation, . Contrast comparison measure maximally if and only if , that is, when the images have equal contrast.

, where is covariance between and , is a structure comparison measure which determines the correlation between the images and . It attains maximal value of if structurally the two images coincide, but its value is equal to zero when there is absolutely no structural coincidence. The quantities , , and are small positive constants that avert the possibility of having zero denominators.

The results of our method were compared to PM method [7], TV method [8], and the D--PM method [11]. Tables 1 and 2, respectively, give a summary of the results from the experiments, using synthetic image in Figure 1 and Lena image in Figure 3. The parameters considered here are thresholding parameter , variance , the time step parameter , and the convolution parameter applied in the D--PM method. The fidelity parameter was dynamically obtained according to (82) or under the AOS scheme in Section 5.1.2, while the rest of the parameters were chosen to guarantee stability and attainment of optimal results.

Table 1: Numerical results for synthetic image .

Algorithm
Parameters
Number of steps
CPU time (sec)
PSNR
MAE
PSNRE
SSIM

PM
30
12
0.25
232
12.73
34.23
2.84
25.26
0.9829

TV
30
n/a
0.2
203
17.86
32.79
3.76
15.76
0.9748

D--PM
30
1
0.25
42
2.40
37.30
3.38
24.76
0.9838

PMS
30
n/a
0.20
262
12.01
35.75
2.15
25.60
0.9870

AOS
30
n/a
3.00
17
2.51
33.30
2.92
25.40
0.9839

Table 2: Numerical results for Lena image .

Algorithm
Parameters
Number of steps
CPU time (sec)
PSNR
MAE
PSNRE
SSIM

PM
30
12
0.25
60
2.43
27.10
7.75
21.70
0.7234

TV
30
n/a
0.2
149
13.87
27.46
7.49
22.63
0.7960

D--PM
30
4
2
13
0.90
28.17
7.35
24.78
0.7880

PMS
30
n/a
0.2
122
5.80
27.73
7.45
26.60
0.8060

AOS
30
n/a
2.0
7
1.24
28.06
7.27
27.80
0.8509

(a) Original image

(b) Noisy image

(c) Proposed model I (PMS)

(d) Proposed model I (AOS)

(e) PM model

(f) TV model

(g) D--PM model

Figure 1: Synthetic image (). (a) Original image. (b) Noisy image corrupted by Gaussian noise for . (c) Our algorithm by PMS; (262 steps). (d) Our algorithm by AOS scheme; (17 steps). (e) PM method; and (232 steps). (f) TV method; (203 steps). (g) D--PM method; , , , and (42 steps).

For nontexture images as displayed in Figure 1 and the results of our method using PMS scheme, as shown in the fourth row of Table 1, demonstrate better performance than those PM and TV as indicated by the higher PSNR, PSNRE, and SSIM values and lower MAE. And, in spite of higher iterative steps, the corresponding CPU time is lower compared to TV and the traditional Perona-Malik (PM) model. And although D--PM model shows better results than PSNR and MAE results, it can be observed that in terms of edge feature recovery measure (PSNRE) and general structural coincidence measure (SSIM) our method performs better. Moreover, implementing our model by the OAS scheme revealed faster execution in terms of both the CPU time and iteration steps (see Table 1). The MAE has lower and even better PSNR results than TV (see Table 1). Looking at visual results of our method using PM scheme (PMS) and AOS, respectively, in Figures 1(c) and 1(d), there is manifestly better visual appeal for our method compared to the PM method (see Figure 3(e)) which shows some speckles, TV method (see Figure 1(f)) which exhibits staircase effects and slight loss in contrast, and the D--method (see Figure 3(g)) which shows slightly deformed edges.

However, for real image such as the given Lena image in Figure 3, the results as shown in Table 2 indicate that our method by AOS scheme gives superior results as shown by the extremely lower iteration steps (7 steps), very short CPU processing time (1.24 sec), better PSNR (28.06), and even lower MAE (7.27) compared to those of the TV method and PM method. Note that our model implemented using AOS performs better than the same model implemented using PM scheme (PMS) for real images. And, in spite of the slightly superior PSNR and MAE values by the D--PM method, our method not only gives better edge preservation as evidenced by the higher PSNRE, but also gives closer structural coincidence than the other three models. Further, the visual appeal of our method, whether using PMS (Figure 3(c)) or AOS scheme (Figure 3(d)), excels those of the traditional Perona-Malik (PM) (see Figure 3(e)) method which shows speckles and a bit of blur, TV method (see Figure 3(f)) which introduces staircasing effects on the denoised image, and even the D--PM method, which shows blockiness and some slight speckle effects.

In addition, from the similarity curves given in Figures 2 and 4, it can be observed in the light of the marked areas, for instance, that this model performs better than its comparisons, both for the synthetic image and real image. Note also that, even though in performance metrics, especially in terms of PSNR and MAE, the D--PM method seems to perform better, the similarity curves attest that our method generates restored images that more closely match the original image than the results obtained by the D--PM method (see Figures 2(e) and 4(e)).

(a) Proposed model I (AOS)

(b) Proposed model I (PMS)

(c) PM model

(d) TV model

(e) D--PM model

Figure 2: Synthetic image (). Similarity graphs between the original image and results of our method (AOS and PMS), PM method, TV method, and D--PM method, respectively.

(a) Original image

(b) Noisy image

(c) Proposed model I (PMS)

(d) Proposed model I (AOS)

(e) PM model

(f) TV model

(g) D--PM model

Figure 3: Lena image (). (a) Original image. (b) Noisy image corrupted by Gaussian noise for . (c) Our algorithm by PMS; (122 steps). (d) Our algorithm by AOS scheme; (7 steps). (e) PM method; and (60 steps). (f) TV method; (149 steps). (g) D--PM method; , , , and (13 steps).

(a) Proposed model I (AOS)

(b) Proposed model I (PMS)

(c) PM model

(d) TV model

(e) D--PM model

Figure 4: Lena image (). Similarity graphs between the original image and results of our method (AOS and PMS), PM method, TV method, and D--PM method, respectively.

6. Conclusion

In this paper, we have proposed a modified total variation model based on the strictly convex modification for image denoising. The main idea was to offer a better image restoration model that is strictly convex and is, therefore, not subject to backward diffusion which has the potential of introducing blurs and to limit the number of parameters available for manual manipulation. This, probably, explains why, in spite of the fact that our model tends to TV as tends to , its practical implementation enjoys better image visual sharpness. This visual sharpness is comparably visible in the results by PM method, but it is still marred behind the oversmoothing and speckle effects in PM images. In fact, the reason why the D--PM method tends to give images that have some blur and deformed edges is attributable to the backward diffusion that it introduces in the course of diffusion. Indeed, it is difficult to assure convergence to a solution of minimum energy given the backward-forward motion of the diffusion process induced by the D--PM algorithm. At the practical level, the increased number of parameters in the D--PM algorithm tends to make it difficult to arrive at a definite permutation of parameter values that would give an optimal result.

For the proposed method, we have demonstrated the existence and uniqueness of the solution of the model. Moreover, numerical implementation of our model also does not suffer potential inaccuracies typical of the TV method, since we do not need to add any small perturbation constant in TV method. Our thresholding parameter depends on the evolution parameter and therefore does not have to be constrained to smaller values as in the case of PM model. The resultant evolution equation has been discretized and implemented using PM scheme and AOS scheme to demonstrate the performance of our algorithm. From the given experimental results, the PSNR values, MAE values, Iterative steps, the CPU processing time, PSNRE, SSIM, and visual appeal of our denoised images all testify that our method is actually a good balance of the PM, TV, and D--PM methods and hence a better image restoration model.

However, in real life application, we observe that the success of any denoising algorithm depends on the type of image being considered, the type of noise, the application intended for the results of the restoration, the extent of degradation, and indeed the implementation platform. And although we have only considered additive Gaussian white noise, different kinds of noise (whether Poisson, speckle, salt, or paper, among others) will require different formulations for the fidelity part and may even demand the application of more than just one formulation for effectiveness. The choice of platform and scheme of implementation must also be appropriately made for efficient performance of the formulation.

With respect to the use of the restoration results, there are situations where the noise removal, generally, may be counterproductive. This usually occurs when the oscillations due to the noise are of comparable scale to those of the features being targeted for preservation. A case like this may require a combination of formulations [42, 43].

For images that are heavily degraded, it might be necessary to do a preconvolution, to blur the noise effects, and then to use an effective formulation such as this one to recover semantically important features such as the edges and contours of the image.

Conflict of Interests

The authors declare that they have no conflict of interests whatsoever and do approve the publication of this paper.

Acknowledgments

This work is partially supported by the National Science Foundation of China (11271100 and 11301113), the Ph.D. Programs Foundation of Ministry of Education of China (no. 20132302120057), the class General Financial Grant from the China Postdoctoral Science Foundation (Grant no. 2012M510933), the Fundamental Research Funds for the Central Universities (Grant nos. HIT. NSRIF. 2011003 and HIT. A. 201412), the Program for Innovation Research of Science in Harbin Institute of Technology (Grant no. PIRS OF HIT A201403), and Harbin Science and Technology Innovative Talents Project of Special Fund (2013RFXYJ044).

References

T. Chang and C. C. J. Kuo, “Texture analysis and classification with tree-structured wavelet transform,” IEEE Transactions on Image Processing, vol. 2, no. 4, pp. 429–441, 1993.

F. Scholkmann, V. Revol, R. Kaufmann, H. Baronowski, and C. Kottler, “A new method for fusion, denoising and enhancement of x-ray images retrieved from Talbot-Lau grating interferometry,” Physics in Medicine and Biology, vol. 59, no. 6, p. 1425, 2014.

J. Ma and G. Plonka, “Combined curvelet shrinkage and nonlinear anisotropic diffusion,” IEEE Transactions on Image Processing, vol. 16, no. 9, pp. 2198–2206, 2007.

E. Candes and D. L. Donoho, “Curvelets: a surprisingly effective nonadaptive representation for objects with edges,” Tech. Rep., 2000, DTIC Document.

E. J. Candès and D. L. Donoho, “New tight frames of curvelets and optimal representations of objects with piecewise singularities,” Communications on Pure and Applied Mathematics, vol. 57, no. 2, pp. 219–266, 2004.

E. Candès, L. Demanet, D. Donoho, and L. Ying, “Fast discrete curvelet transforms,” Multiscale Modeling & Simulation, vol. 5, no. 3, pp. 861–899, 2006.

P. Perona and J. Malik, “Scale-space and edge detection using anisotropic diffusion,” IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 12, no. 7, pp. 629–639, 1990.

L. I. Rudin, S. Osher, and E. Fatemi, “Nonlinear total variation based noise removal algorithms,” Physica D, vol. 60, no. 1–4, pp. 259–268, 1992.

T. F. Chan and S. Esedo\=glu, “Aspects of total variation regularized function approximation,” SIAM Journal on Applied Mathematics, vol. 65, no. 5, pp. 1817–1837, 2005.

J. Weickert, B. M. Ter Haar Romeny, and M. A. Viergever, “Efficient and reliable schemes for nonlinear diffusion filtering,” IEEE Transactions on Image Processing, vol. 7, no. 3, pp. 398–410, 1998.

Z. Guo, J. Sun, D. Zhang, and B. Wu, “Adaptive Perona-Malik model based on the variable exponent for image denoising,” IEEE Transactions on Image Processing, vol. 21, no. 3, pp. 958–967, 2012.

F. Catté, P.-L. Lions, J.-M. Morel, and T. Coll, “Image selective smoothing and edge detection by nonlinear diffusion,” SIAM Journal on Numerical Analysis, vol. 29, no. 1, pp. 182–193, 1992.

Y. Chen and M. Rao, “Minimization problems and associated flows related to weighted energy and total variation,” SIAM Journal on Mathematical Analysis, vol. 34, no. 5, pp. 1084–1104, 2003.

B. J. Maiseli, Q. Liu, O. A. Elisha, and H. Gao, “Adaptive Charbonnier superresolution method with robust edge preservation capabilities,” Journal of Electronic Imaging, vol. 22, no. 4, Article ID 043027, 2013.

J. Shah, “Common framework for curve evolution, segmentation and anisotropic diffusion,” in Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition (CVPR '96), pp. 136–142, June 1996.

D. Barash and D. Comaniciu, “A common framework for nonlinear diffusion, adaptive smoothing, bilateral filtering and mean shift,” Image and Vision Computing, vol. 22, no. 1, pp. 73–81, 2004.

J. J. Koenderink, “The structure of images,” Biological Cybernetics, vol. 50, no. 5, pp. 363–370, 1984.

A. P. Witkin, “Scale-space filtering,” US patent no. 4,658,372, 1987

J. Weickert, Anisotropic Diffusion in Image Processing, vol. 1, Teubner, Stuttgart, Germany, 1998.

Y. Chen and T. Wunderli, “Adaptive total variation for image restoration in BV space,” Journal of Mathematical Analysis and Applications, vol. 272, no. 1, pp. 117–137, 2002.

T. F. Chan and J. Shen, “Mathematical models for local nontexture inpaintings,” SIAM Journal on Applied Mathematics, vol. 62, no. 3, pp. 1019–1043, 2001.

C. R. Vogel, “Total variation regularization for Ill-posed problems,” Tech. Rep., Department of Mathematical Sciences, Montana State University, 1993.

L. Vese, Problemes variationnels et EDP pour lA analyse dA images et lA evolution de courbes [Ph.D. thesis], Universite de Nice Sophia-Antipolis, 1996.

Y. Chen, S. Levine, and M. Rao, “Variable exponent, linear growth functionals in image restoration,” SIAM Journal on Applied Mathematics, vol. 66, no. 4, pp. 1383–1406, 2006.

T. Chan, A. Marquina, and P. Mulet, “High-order total variation-based image restoration,” SIAM Journal on Scientific Computing, vol. 22, no. 2, pp. 503–516, 2000.

F. Andreu-Vaillo, V. Caselles, and J. M. Mazón, Parabolic QuasiLinear Equations Minimizing Linear Growth Functionals, vol. 223, Springer, 2004.

R. Acar and C. R. Vogel, “Analysis of bounded variation penalty methods for ill-posed problems,” Inverse Problems, vol. 10, no. 6, pp. 1217–1229, 1994.

D. M. Strong and T. F. Chan, “Spatially and scale adaptive total variation based regularization and anisotropic diffusion in image processing,” in Diusion in Image Processing, UCLA Math Department CAM Report, Cite-seer, 1996.

A. Chambolle and P.-L. Lions, “Image recovery via total variation minimization and related problems,” Numerische Mathematik, vol. 76, no. 2, pp. 167–188, 1997.

Y.-L. You, W. Xu, A. Tannenbaum, and M. Kaveh, “Behavioral analysis of anisotropic diffusion in image processing,” IEEE Transactions on Image Processing, vol. 5, no. 11, pp. 1539–1553, 1996.

L. Vese, “A study in the BV space of a denoising-deblurring variational problem,” Applied Mathematics and Optimization, vol. 44, no. 2, pp. 131–161, 2001.

A. Marquina and S. Osher, “Explicit algorithms for a new time dependent model based on level set motion for nonlinear deblurring and noise removal,” SIAM Journal on Scientific Computing, vol. 22, no. 2, pp. 387–405, 2000.

X. Zhou, “An evolution problem for plastic antiplanar shear,” Applied Mathematics and Optimization, vol. 25, no. 3, pp. 263–285, 1992.

G. Aubert and P. Kornprobst, Mathematical Problems in Image Processing, vol. 147, Springer, New York, NY, USA, 2nd edition, 2006.

J. L. Doob, Measure Theory, vol. 143, Springer, New York, NY, USA, 1994.

L. Ambrosio, N. Fusco, and D. Pallara, Functions of Bounded Variation and Free Discontinuity Problems, vol. 254, Clarendon Press, Oxford, UK, 2000.

C. P. Niculescu and L.-E. Persson, Convex Functions and Their Applications: : A Contemporary Approach, vol. 23, Springer, 2006.

M. Renardy and R. C. Rogers, An Introduction to Partial Differential Equations, vol. 13, Springer, 2nd edition, 2004.

H. Brézis, Opérateurs Maximaux Monotones et Semi-Groupes de Contractions dans les Espaces de Hilbert, vol. 50, North-Holland, 1973.

S. Durand, J. Fadili, and M. Nikolova, “Multiplicative noise removal using fidelity on frame coefficients,” Journal of Mathematical Imaging and Vision, vol. 36, no. 3, pp. 201–226, 2010.

Z. Wang, A. C. Bovik, H. R. Sheikh, and E. P. Simoncelli, “Image quality assessment: from error visibility to structural similarity,” IEEE Transactions on Image Processing, vol. 13, no. 4, pp. 600–612, 2004.

E. A. Ogada, Z. Guo, and B. Wu, “An alternative variational framework for image denoising,” Abstract and Applied Analysis, vol. 2014, Article ID 939131, 16 pages, 2014.

Y. Yu and S. T. Acton, “Speckle reducing anisotropic diffusion,” IEEE Transactions on Image Processing, vol. 11, no. 11, pp. 1260–1270, 2002.