Abstract and Applied Analysis

Volume 2013 (2013), Article ID 585310, 10 pages

http://dx.doi.org/10.1155/2013/585310

Research Article

Fractional-Order Total Variation Image Restoration Based on Primal-Dual Algorithm

Dali Chen,¹ YangQuan Chen,² and Dingyu Xue¹

¹College of Information Science and Engineering, Northeastern University, Shenyang, Liaoning 110006, China

²MESA Lab, University of California, Merced, 5200 North Lake Road, Merced, CA 95343, USA

Received 7 August 2013; Accepted 27 September 2013

Academic Editor: Dumitru Baleanu

Copyright © 2013 Dali Chen 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.

This paper proposes a fractional-order total variation image denoising algorithm based on the primal-dual method, which provides a much more elegant and effective way of treating problems of the algorithm implementation, ill-posed inverse, convergence rate, and blocky effect. The fractional-order total variation model is introduced by generalizing the first-order model, and the corresponding saddle-point and dual formulation are constructed in theory. In order to guarantee convergence rate, the primal-dual algorithm was used to solve the constructed saddle-point problem, and the final numerical procedure is given for image denoising. Finally, the experimental results demonstrate that the proposed methodology avoids the blocky effect, achieves state-of-the-art performance, and guarantees convergence rate.

1. Introduction

Since the work of Rudin et al. [1], total variation (TV) minimization problems arise in many image processing applications for regularizing inverse problems where one expects the recovered image or signal to be piecewise constant [2–5]. The typical total variation model, for example, ROF model, has been proved to be able to achieve a good tradeoff between edge preservation and noise removal [1]. However, it tends to produce the so-called blocky (staircase) effects on the images because it favors a piecewise constant solution in bounded variation (BV) space [6]. In order to deal with blocky effects, the modification of TV model, which generalizes the differential order in regularization term, has aroused the more and more attentions of numerous scholars.

The improved methods of TV model are divided into two kinds: the high-order derivative and the fractional-order derivative. The first one replaces the first-order derivative in regularization term by the high-order derivative. For example, a fourth-order partial differential equation-(PDE-based) denoising mode was proposed by [7], in which the regularized solution is obtained by solving the minimization of potential function of second-order derivative of the image. It has been proved that this model is able to deal with the blocky problem. However, it tends to cause the sign of uplifting effect and formation of artifacts around edges. For this problem, an improved fourth-order PDE model, replacing the Laplacian operator of diffusivity function by the gradient operator, was proposed in [8].

In this paper, our interest focuses on the second generalization which deals with fractional-order differentiation. Fractional calculus is a rapidly growing mathematical discipline, which provides an important tool for nonlocal field theories [9]. Recently, it has been greatly studied in computer vision [10–16], and the main reason for this development is the expectation that the use of this theory will lead to a much more elegant and effective way of treating problems of blocky effect and detailed information protection. Specially, the fractional-order total variation (TV) models play an important role for image denoising, inpainting, and motion estimation [17–21]. So far, the methods adopted to deal with the fractional-order total variation problem are divided into three kinds: solving the associated Euler-Lagrange equation, which is a nonlinear partial differential equation (PDE) [17, 18]; using the methods based on duality [6]; using the method based on majorization-minimization (MM) algorithm [22]. Although the convergence of these methods has been studied in the corresponding literatures, the rate of convergence for them has not been discussed yet. For this problem, in this paper, a primal-dual algorithm is used to solve the fractional-order total variation problem, which is able to guarantee convergence rate.

The primal-dual algorithm was first presented by [23] and named as the primal-dual hybrid gradient (PDHG) method in [24]. In this method, each iteration updates both a primal and a dual variable. It is thus able to avoid some of the difficulties that arise when working only on the primal or dual side [25, 26]. The convergence of the primal-dual algorithm has been studied in [27, 28]. Recently a unified form of primal-dual algorithm was presented by [29, 30], which demonstrated that, in some case, these algorithms can also achieve the rate of convergence. In our paper, a new image denoising method is proposed, in which the primal-dual algorithm is used to solve the fractional-order total variation denoising model. The proposed denoising method is able to avoid the blocky effect, achieves state-of-the-art performance, and guarantees the rate of convergence.

This paper is organized as follows. Section 1 introduces prior work, focusing on the main problems with existing methods that are addressed by our model. In Section 2, the fractional-order total variation denoising model is described, and the corresponding saddle-point and dual formulation are constructed in theory. Based on this, the primal-dual algorithm was used to solve the constructed saddle-point problem, and the final numerical implementation is presented for image denoising. Experimental evaluation is presented in Section 3 and the paper is concluded in Section 4.

2. Fractional-Order Total Variation Model and Primal-Dual Algorithm

2.1. Model Description

Let denote the observed noisy image, where denotes the location with a rectangular image domain and is white Gaussian noise. The typical total variation (TV) denoising model estimates the desired clean image by solving the following finite-dimensional optimization problem:

where is the regularization term, is the data fidelity term, is -norm of , and is regularization parameter which controls the degree of smoothing. In this paper, we consider a fractional-order total variation model, defined as

where is obtained by the following formula:

where is the discretization of , is the fractional-order derivative operator, and . and are linear operators corresponding to horizontal and vertical fractional-order derivative.

From Grünwald-Letnikov fractional derivative definition [31], the finite fractional-order forward difference can be obtained by

where , denotes the generalized binomial coefficient, and is the Gamma function. In addition, the coefficients can also be obtained recursively from

When , for and (4) is the first-order forward derivative as usual.

To describe (2) in matrix algebra language, we reord the image matrix and row wisely into the vector and , associating the element of the two-dimensional structure with the element of the vector structure, , and . We have and , where . The component of the fractional-order derivative can thus be represented as a multiplication of the vector by a matrix , for :

From this definition, the discrete version of the prime fractional-order total variation (FOTV) model (2) can be written as

Define a norm where , and (7) can be rewritten in term of as

Here .

2.2. Dual Problem and Saddle-Point Problem

For any vector, and . We have , where is the dual norm of norm , and the and are analogous to and , respectively.

Let the function be the conjugate of ; we have

in which case we can interpret as the indicator function for the unit ball in the dual norm, and with . Since if is closed and convex, we have

Substituting (10) into (8), we can obtain the saddle-point formulation of FOTV model defined by

The minimization problem of (11) can be solved exactly as

Substituting (12) into (11) yields the following dual problem:

where is the corresponding discretization of the fractional-order divergence, and when , it is the discrete divergence as usual.

For primal-dual feasible pair , the partial primal-dual gap is defined by

If is the saddle point of problem (11), then is primal optimal, is dual optimal [32], and we have

which vanishes only if is itself the saddle point [29]. The primal-dual gap is a measure of closeness of the primal-dual to the primal-dual solution, and we use it to design the stopping criterion for our numerical algorithm in this paper.

2.3. Primal-Dual Method for Solving Saddle-Point Problem

The first-order primal-dual method summarized in [29] for convex problem was applied in this paper to solve the saddle-point problem described by (11), which is described by Algorithm 1.

Step 1. Choose with , , , and .

Step 2. Update , , , , , as follows:

Algorithm 1

Let and ; then there exists a saddle-point such that and . Since is uniformly convex such that it has a Lipschitz continuous gradient, the output possesses a rate of convergence of [29]. More specifically, we assume that the constant stepsizes are used; that is, , , and for all . If ; then the convergence rate can be guaranteed for , and convergence rate can be guaranteed for [30].

To summarize, our entire noise removal algorithm in a form of a pseudocode is done in following:

(1) initialization:

(i)for a given noisy image of size , we reord the image row wisely into the vector with , ; (ii)compute the fractional-order discrete operator according to (6); (iii)initialize the regularization parameter , iteration number , and maximum permissible error ; (iv)set , , , , , and ;

(2) iteration: compute by the following steps.

Step 1. Solve the problem (A) by

Step 2. Substitute into the (B) and solve the problem (B) by

Step 3. Compute , , and .

Step 4. Update by .

Step 5. Compute the primal-dual gap according to (15) as follows:

Step 6. If () or (), then we terminate the iteration and output ; otherwise, go back to Step 1.

Due to the simplicity of our model, it is easy to be implemented and can be effectively accelerated on parallel hardware such as field programmable gata array (FPGA) and graphics processing unit (GPU). In order to employ the huge computational power and the parallel processing capabilities of FPGA to obtain a fully accelerated implementation of our denoising method, the numerical method should be working on the regular grids. For this purpose, the implementation of our algorithm on FPGA is described as the following strategy. First, an image frame is captured and read in the processor from a live camera. For the captured frame , the th iteration of our method performs the following updates on prime variable and dual variable :

where the dual variable corresponds with and the discrete fractional divergence operator is adjoint to the discrete fractional gradient operator . Note that (19) and (20), respectively, correspond with the solutions of problem (A) and (B). The updates of the other parameters are the same as the description in the foregoing pseudocode, and the number of iteration is fixed to ensure the correct timing sequence. After the iteration, the denoised image is outputted, the next image frame is captured, and the same process is repeated until the assignment is over.

3. Experiments and Analysis

3.1. Restraint of Block Effect

The blocky effect is the main drawback produced by the typical first-order total variation denoising algorithm. In this section, some experiments are given to assess the capability of reducing block effect of our proposed fractional-order total variation (FOTV) denoising algorithm. Firstly, the one-dimensional signal is used as the test signal, . We contaminate the given signal using the additive white Gaussian noise (AWGN) with standard deviation (SD) of . The corrupted signal is shown in Figure 1(a). The proposed FOTV denoising algorithm with and is used to process the contaminated signal, respectively. When , our FOTV algorithm is the typical ROF denoising algorithm based on primal-dual method. In this experiment, we set , and the results at the 1000th iteration are, respectively, shown in Figures 1(b) and 1(c). The difference is obvious: while the first-order TV denoising algorithm approximates the observed signal with a step signal, the fractional-order algorithm with a piecewise planar signal which looks more natural and does not produce false edges.

(a) Noisy signal

(b) Processed by FOTV with

(c) Processed by FOTV with

Figure 1: Signal processed by FOTV with and .

We now consider a two-dimensional image “Lena” of size . The image is degraded by using AWGN with standard deviation of and the result is shown in Figure 2(a). This degraded image is, respectively, fed into the proposed FOTV denoising algorithm with and as initial condition, and the time evolution of the algorithms begins. The results at the 1000th iteration are, respectively, shown in Figures 2(b) and 2(c). The blocky effects are obvious in Figure 2(b), while the Figure 2(c) looks more natural and does not produce blocky effect.

(a) Noisy image

(b) Processed by FOTV with

(c) Processed by FOTV with

Figure 2: Image processed by FOTV with and .

In conclusion, the fractional-order TV algorithm can reduce blocky effect effectively comparing with the traditional first-order TV algorithm.

3.2. Analysis of Denoising Performance

The aim of this section is to analyze the denoising performance of the FOTV algorithm. For this purpose, three famous test images, called Barbara, Lena, and Peppers, are used in the experiments, as shown in Figure 3. In order to quantify the denoised image, we consider the peak-signal-to-noise ratio (PSNR), which has been largely used in the literature and commonly applied to determine the quality of a processed image. It can be calculated by the following formula:

where is the original image, is the denoised image, and is the maximal gray level of the image.

(a) Barbara

(b) Lena

(c) Peppers

Figure 3: Original image: Barbara, Lena, and Peppers.

Firstly, in order to decide the value of the fractional-order , we study the relation between the PSNR and . Figure 4 shows the relation between the PSNR and on Lena and Peppers image corrupted by the additive white Gaussian noise with standard deviation (SD) of . The left image is the result of Lena image and the right one is the result of Peppers image. From the figure, we can obtain the following conclusions. Firstly, the PSNR reaches a maximum between and . Secondly, the PSNR at is lower than PSNR at , which is owing to the blocky effect. Finally, the PSNR decreases rapidly as tends to zero, and when , the PSNR reaches a minimum. According to these facts, we can select the fractional-order between and .

(a) Lena

(b) Peppers

Figure 4: Relation between the fractional-order and PSNR on Lena (a) and Peppers (b) image corrupted by the additive white Gaussian noise with standard deviation (SD) of .

In order to analyze the denoising performance, the additive white Gaussian noise (AWGN) with standard deviation (SD) of , , and is, respectively, added in the three test images. Five denoising algorithms are used to process these given noisy images, which are the improved Perona and Malik (IP-M) model [33], fourth order (F-O) PDE model [7], improved fourth order (IF-O) PDE model [8], ROF model [1], and the proposed FOTV denoising algorithm. In this experiment, we set and for IP-M, F-O-PDE, and IF-O-PDE models, and we set , , and for ROF model. These parameters are able to ensure the best denoising performance of the corresponding denoising model. Based on the conclusion of the foregoing experiment, we set for FOTV denoising algorithm, and the denoising results are shown in Table 1. In the table, the first column lists the given image and the second column lists the standard deviation of noisy image. The PSNRs of the denoised images processed by the different denoising algorithms are listed under the corresponding denoising model. The bigger the PSNR is, the better the denoising performance is. For easy observation, the biggest PSNR values are shown in boldface. From the table, it is obvious that the PSNR of our proposed algorithm is bigger than that of the other four algorithms, so we can conclude that our FOTV denoising algorithm outperforms the other algorithms.

Table 1: PSNR quantitative comparison among five denoising models.

Image
SD
IP-M
F-O-PDE
IF-O-PDE
ROF
FOTV

Barbara
10
31.2427
29.3776
29.3777
31.0871
31.2850

20
26.6231
24.8155
24.8156
26.8212
27.0976

30
24.4496
22.9267
22.9273
24.7415
24.9479

Lena
10
33.6393
31.5440
31.5442
33.8378
34.4747

20
29.7665
27.8080
27.8098
30.4101
31.2153

30
27.4437
26.0883
26.1019
28.5889
29.3855

Peppers
10
33.6967
31.8199
31.8205
33.8715
34.0619

20
30.0275
28.0437
28.0457
30.1768
31.0285

30
27.4982
26.1613
26.1694
28.2689
29.0046

In order to further verify the denoising performance of our proposed algorithm, Figure 5 shows the denoised images processed by four denoising algorithms. The first figure is the famous image “Barbara” corrupted by the additive white Gaussian noise with standard deviation (SD) of . In order to show the comparison clearly, the partial enlarged view of the noisy image is shown in the second figure. The denoised result of IP-M algorithm is shown in the right-top figure. The second row shows the denoised images processed by IF-O PDE algorithm, ROF algorithm, and our FOTV algorithm. In Figures 5(c) and 5(d), although the detailed information is preserved, there are a lot of noises unremoved. The result in Figure 5(e) looks blocky. Only the result of our model, as shown in Figure 5(f), looks natural and does not produce false edges. So we conclude that our proposed FOTV algorithm is able to achieve a better tradeoff between edge preservation and noise removal.

(a) Barbara with 20 SD Gaussian noise

(b) Partial view of (a)

(c) IP-M algorithm

(d) IF-O PDE algorithm

(e) ROF algorithm

(f) FOTV algorithm

Figure 5: Comparison of four algorithms for the noisy image with additive Gaussian white noise of .

3.3. The Convergence and the Rate of Convergence

The convergence and the rate of convergence are two important factors for evaluating the performance of the denoising method. In this section, some experiments are given to show the convergence and convergence rate of our proposed fractional-order total variation (FOTV) denoising algorithm.

For this purpose, the proposed FOTV algorithm is, respectively, used to process the image “Lena” corrupted by the additive white Gaussian noise with standard deviation (SD) of 10, 20, and 30, and the PSNR is recorded at each iteration. Figure 6(a) shows the relation between the iteration and PSNR. The horizontal axis of this figure is iteration number and the vertical axis is PSNR. The different colorful curve denotes the result obtained from the different noisy image. The red dotted line is the result of the noisy image with , the green dashed line is the result of the noisy image with , and the red real line is the result of the noisy image with . From the figure, we can see that the PSNR values are stable after 50 iterations, which indicate that our proposed FOTV denoising algorithm is convergent. In addition, it can be seen that the bigger the variance of white noise is, the longer the stopping times are and the lower the PSNRs are.

(a) Relation between the iteration and PSNR

(b) Relation between the iteration and primal-dual gap

Figure 6: Relation among the iteration, PSNR, and primal-dual gap on the noisy Lena image with , , and .

Figure 6 plots the convergence of the FOTV algorithm for the above experiment together with the theoretical and rate. For easy observation, the horizontal axis is the log curve of iteration number and the vertical axis is the log curve of primal-dual gap calculated in Step 5. The dashed line shows the theoretical rate and the dotted line shows the theoretical rate. The red line shows the convergence of FOTV algorithm in the noisy image with , the green line shows the convergence of FOTV algorithm in the noisy image with , and the blue line shows the convergence of FOTV algorithm in the noisy image with . From this figure, it can be seen that the smaller the variance of white noise is, the faster convergence is. In addition, we can see that the FOTV algorithm converges at rate firstly and rate in the end.

4. Conclusion

A fractional-order total variation image denoising algorithm based on primal-dual method was proposed in this paper. The main contributions are as follows: the saddle point and dual formulation of the proposed fractional-order total variation model are constructed in theory; the primal-dual algorithm was used to solve the constructed saddle-point problem, in which convergence rate can be guaranteed; the experiments demonstrate that the proposed methodology is able to avoid the blocky effect, achieve state-of-the-art performance, and guarantee convergence rate.

Acknowledgments

This work was supported by National Natural Science Foundation of China (no. 61201378), Scientific Research Fund of Liaoning Provincial Education Department (L2012073), and Fundamental Research Funds for the Central Universities (N110304001).

References

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

J.-F. Aujol, “Some first-order algorithms for total variation based image restoration,” Journal of Mathematical Imaging and Vision, vol. 34, no. 3, pp. 307–327, 2009.

C. Vogel and M. Oman, “Fast, robust total variation-based reconstruction of noisy, blurred images,” IEEE Transactions on Image Processing, vol. 7, no. 6, pp. 813–824, 1998.

F. Alter, S. Durand, and J. Froment, “Adapted total variation for artifact free decompression of JPEG images,” Journal of Mathematical Imaging and Vision, vol. 23, no. 2, pp. 199–211, 2005.

F. Li, C. Shen, and C. Li, “Multiphase soft segmentation with total variation and H ¹ regularization,” Journal of Mathematical Imaging and Vision, vol. 37, no. 2, pp. 98–111, 2010.

J. Zhang, Z. Wei, and L. Xiao, “Adaptive fractional-order multi-scale method for image denoising,” Journal of Mathematical Imaging and Vision, vol. 43, no. 1, pp. 39–49, 2012.

Y. L. You and M. Kaveh, “Fourth-order partial differential equations for noise removal,” IEEE Transactions on Image Processing, vol. 9, no. 10, pp. 1723–1730, 2000.

M. Hajiaboli, “A self-governing fourth-order nonlinear diffusion filter for image noise removal,” IPSJ Transactions on Computer Vision and Applications, vol. 2, pp. 94–103, 2010.

R. Herrmann, Fractional Calculus: An Introduction for Physicists, World Scientific, New Jersey, NJ, USA, 2011.

S.-C. Liu and S. Chang, “Dimension estimation of discrete-time fractional Brownian motion with applications to image texture classification,” IEEE Transactions on Image Processing, vol. 6, no. 8, pp. 1176–1184, 1997.

D. Chen, Y. Chen, and D. Xue, “Digital fractional order Savitzky-Golay differentiator,” IEEE Transactions on Circuits and Systems II, vol. 58, no. 11, pp. 758–762, 2011.

S. Didas, B. Burgeth, A. Imiya, and J. Weickert, “Regularity and scale-space properties of fractional high order linear filtering,” in 5th International Conference on Scale Space and PDE Methods in Computer Vision, Scale-Space 2005, pp. 13–25, deu, April 2005.

B. Ninness, “Estimation of Noise,” IEEE Transactions on Information Theory, vol. 44, no. 1, pp. 32–46, 1998.

I. Petras, D. Sierociuk, and I. Podlubny, “Identification of parameters of a half-order system,” IEEE Transactions on Signal Processing, vol. 60, no. 10, pp. 5561–5566, 2012.

D. Chen, Y. Chen, and D. Xue, “1-D and 2-D digital fractional-order Savitzky-Golay differentiator,” Signal, Image and Video Processing, vol. 6, no. 3, pp. 503–511, 2012.

Y. F. Pu, J. L. Zhou, and X. Yuan, “Fractional differential mask: a fractional differential-based approach for multiscale texture enhancement,” IEEE Transactions on Image Processing, vol. 19, no. 2, pp. 491–511, 2010.

D. Chen, H. Sheng, Y. Chen, and D. Xue, “Fractional-order variational optical flow model for motion estimation,” Philosophical Transactions of the Royal Society A, vol. 371, no. 1990, 2013.

J. Bai and X. Feng, “Fractional-order anisotropic diffusion for image denoising,” IEEE Transactions on Image Processing, vol. 16, no. 10, pp. 2492–2502, 2007.

P. Guidotti and J. V. Lambers, “Two new nonlinear nonlocal diffusions for noise reduction,” Journal of Mathematical Imaging and Vision, vol. 33, no. 1, pp. 25–37, 2009.

E. Cuesta, M. Kirane, and S. A. Malik, “Image structure preserving denoising using generalized fractional time integrals,” Signal Processing, vol. 92, no. 2, pp. 553–563, 2012.

M. Janev, S. Pilipović, T. Atanacković, R. Obradović, and N. Ralević, “Fully fractional anisotropic diffusion for image denoising,” Mathematical and Computer Modelling, vol. 54, no. 1-2, pp. 729–741, 2011.

D. Chen, S. Sun, C. Zhang, Y. Chen, and D. Xue, “Fractional-order TV-L ² model for image denoising,” Central European Journal of Physics, 2013.

K. Arrow, L. Hurwicz, and H. Uzawa, “Studies in linear and non-linear programming,” in Stanford Mathematical Studies in the Social Sciences, Stanford University Press, Palo Alto, Calif, USA, 1958.

M. Zhu and T. Chan, “An efficient primal-dual hybrid gradient algorithm for total variation image restoration,” CAM Report 08-34, UCLA, Los Angeles, Calif, USA, 2008.

A. Chambolle, “An Algorithm for Total Variation Minimization and Applications,” Journal of Mathematical Imaging and Vision, vol. 20, no. 1-2, pp. 89–97, 2004.

E. Esser, X. Zhang, and T. F. Chan, “A general framework for a class of first order primal-dual algorithms for convex optimization in imaging science,” SIAM Journal on Imaging Sciences, vol. 3, no. 4, pp. 1015–1046, 2010.

S. Bonettini and V. Ruggiero, “On the convergence of primal-dual hybrid gradient algorithms for total variation image restoration,” Journal of Mathematical Imaging and Vision, vol. 44, no. 3, pp. 236–253, 2012.

B. He and X. Yuan, “Convergence analysis of primal-dual algorithms for a saddle-point problem: from contraction perspective,” SIAM Journal on Imaging Sciences, vol. 5, no. 1, pp. 119–149, 2012.

A. Chambolle and T. Pock, “A first-order primal-dual algorithm for convex problems with applications to imaging,” Journal of Mathematical Imaging and Vision, vol. 40, no. 1, pp. 120–145, 2011.

Y. Chen, G. Lan, and Y. Ouyang, “Optimal primal-dual methods for a class of saddle point problems,” CAM Report 13-31, UCLA, Los Angeles, Calif, USA, 2013.

I. Podlubny, Fractional Differential Equations, Academic Press, New York, NY, USA, 1999.

S. Boyd and L. Vandenberghe, Convex Optimization, Cambridge University Press, New York, NY, USA, 2004.

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.