Abstract
We address the problem of speckle noise removal. The classical total variation is extensively used in this field to solve such problem, but this method suffers from the staircase-like artifacts and the loss of image details. In order to resolve these problems, a nonconvex total generalized variation (TGV) regularization is used to preserve both edges and details of the images. The TGV regularization which is able to remove the staircase effect has strong theoretical guarantee by means of its high order smooth feature. Our method combines the merits of both the TGV method and the nonconvex variational method and avoids their main drawbacks. Furthermore, we develop an efficient algorithm for solving the nonconvex TGV-based optimization problem. We experimentally demonstrate the excellent performance of the technique, both visually and quantitatively.
1. Introduction
Images generated by coherent imaging modalities, for example, synthetic aperture radar (SAR), ultrasound, and laser imaging, inevitably come with multiplicative noise (also known as speckle), due to the coherent nature of the scattering phenomena. The speckle noise seriously interferes with the upper tasks, such as object recognition [1] and image segmentation [2]. Due to the coherent nature of the image acquisition process, in the speckle noise models, the noise field is multiplied by (not added to) the original image, and it is described by a non-Gaussian probability density function, with Rayleigh and Gamma being common models [3]. So it is signal independent, non-Gaussian, and spatially dependent. Hence, speckle denoising is a very challenging problem compared with additive Gaussian noise.
Speckle noise removal methods have been discussed in many references. Popular methods include bilateral filtering for despeckling [4], wavelet based despeckling approaches [5], and nonlocal means (NL-means) [6]. We will focus on the variational approach for speckle noise removal.
To the best of our knowledge, there exist several variational approaches devoted to speckle noise removal problem, which minimize some appropriate energy functionals, composing a regularization term and a data fitting term. The first total variation-based speckle noise removal method (RLO-method) was presented by Osher et al. [7], which used a constrained optimization approach with two Lagrange multipliers. Aubert and Aujol [8] propose their speckle noise removal method (AA method) in the framework of the maximum a posteriori probability (MAP) estimation. In [9], Huang et al. proposed a new total variation (TV) method for speckle noise removal based on the Aubert-Aujol (AA) method. Because of the nonconvexity of AA method, the global solution is hard to find. To resolve this problem, Bioucas-Dias and Figueiredo [10] apply the MAP estimation method in the log domain and propose a convex speckle noise removal method (BF method). In addition, Steidl and Teuber [11] also propose a convex method (ST method), in which the I-divergence is used as the fidelity term; the reason for this is that the Euler-Lagrange equation of the ST method is equivalent to that of the BF method in the sense of log transform [11].
The methods mentioned above have a common feature that is using the convex TV regularization, which yields piecewise smooth estimates adapted to the structure of the underlying reflectance. Indeed, solutions of variational problems with TV regularization admit many desirable properties, most notably the appearance of sharp edges. However, the regularization with TV also has the so-called staircasing artifacts in the smooth image regions. To overcome the drawback, the total generalized variation (TGV) regularization [12] also has been investigated in the recent work [13], which incorporates the TGV penalty into the existing data fidelity term for the speckle noise removal, and develops two novel variational despeckling methods. TGV-based despeckling method outperforms the traditional TV methods by reducing the staircasing artifacts.
Recently, the developments of nonconvex world involving a variety of applications show that nonconvex regularizer has advantages over convex regularization for restoring images with neat edges. Nikolova et al. [14] discuss the properties of the nonconvex regularizer for removing additive noise and then [15, 16] study data fitting and concave regularization for image recovery problem. Han et al. [17] apply nonconvex TV regularizer to the speckle noise removal, which better preserves edges of restored images compared to classical TV regularizer-based methods. Ochs et al. [18] combine nonconvexity with total generalized variation for reducing the additive noise, deconvolution, and some other applications.
In this study, we incorporate the nonconvex TGV penalty into the existing data fidelity term for speckle removal and develop a novel variation despeckling method. A nonconvex TGV regularizer is used to make the whole regions of the image efficiently smoothed while it is used to make the edges well preserved. As demonstrated in our numerical experiments, the nonconvex TGV-based despeckling method not only outperforms the TGV methods (TGVSNR algorithm) [13] by better preserving edges of images but also is far better than the nonconvex TV-based method (NRSNR algorithm) [17] by removing the staircasing artifacts.
The rest of this paper is organized as follows. In Section 2, some related works are reviewed and discussed. In Section 3, we present a brief review of nonconvex TGV. Then, the new variational method for speckle removal is proposed. In Section 4, a fast algorithm corresponding to the new method is designed. Subsequently, Section 5 describes the experiment results. In Section 6, the final conclusions are drawn.
2. Related Work and Discussion
Let denote an -pixel observed image, assumed to be a sample of a random image . It is known that can be assumed to be the product of the underlying true image intensity and the speckle noise : The probability density function of for the -look SAR image is given by the following Gamma distribution: where is the usual Gamma function. Several major TV-based variational methods are presented below.
The MAP criterion is applied to (1); Aubert and Aujol [8] derived a new method (AA method): where denotes the despeckled result and denotes the regularization parameter.
From the minimization problem (3), we can observe that the global solution of the AA method is hard to find because it has a nonconvex fidelity term. To resolve this problem, the authors in [10] take logarithmic transformation to convert (1) into an additive form:
For simplicity, the pixel subscript of (4) has been dropped. The probability density of the random variable is
Then, according to conditional independence assumption [10], the following equation is derived: where is a constant. Using the MAP criterion, the restored image of the BF method can be inferred by solving the following variational problem:
However, since the logarithmic transformation is nonlinear, which results in luminance distortion of the image, to avoid this undesirable defect, Steidl and Teuber [11] showed that the classical I-divergence method given in (8) does not require nonlinear transformation and shares the same solution of method (7) theoretically:
3. Problem Formulation
In this section, we begin with a brief tour of nonconvex TGV and then propose a new variational method based on it.
3.1. Nonconvex Total Generalized Variation
First, recall the definition of the total generalized variation. For convenience, we assume that is a nonempty, open, and connected set; here, , , is a fixed space dimension. In this paper, we have for the images. The TGV of order [12] and positive weights is defined as where denotes the space of symmetric tensors of order with arguments in . denotes the vector space of compactly supported symmetric tensor field. Note that is a seminorm that equals 0 for all polynomials of degree less than . Consequently, the image restoration with TGV regularization brings about piecewise polynomial intensities.
Specifically, we use the second-order TGV based on the above definition throughout the paper:
Definition (10) is generalized to represent minimization problem itself via Fenchel duality theory [12]:
Here, denotes the symmetric derivative. Such a definition provides a way of balancing between the first and second derivatives of a function. More details on TGV can be found in [12, 19].
Nonconvex variant of TGV, that is, nonconvex TGV, is introduced [18] by Ochs et al. who put forward this theory, which combines with total generalized variation and nonconvex regulation: that is, the regularized -norm, . The nonconvex regulation based upon the TGV can retain all the advantages of both the TGV and nonconvex property. In particular, the despeckling with nonconvex TGV regularization both can lead to piecewise polynomial smoothness and can preserve neat edges of the images. Experiment results in Section 5 can illustrate these facts.
3.2. New Nonconvex Variational Method
Following the MAP estimation process as BF method [10], we propose a novel variational method based on for removing speckle noise:
Let the sum of data term and nonconvex TGV regularization term in formula (13) denote the following:
The regularization term of formula (13) has the advantage of protecting edges of images, but because of the nonconvexity and nonlinearity of the energy function in method (13) classical optimization algorithms cannot be directly used to solve this problem. In the next section, we will design a new fast iteration algorithm.
4. Explicit Numerical Scheme
In the following, to compute the minimizer of (13), we first convert the unconstrained problem (13) into an equivalent constrained one in order to resolve the equivalent problem into two subproblems by the fast augmented Lagrange multiplier (ALM) method. For the first subproblem, as its convexity, the Newton method can be used to solve it. As the second subproblem is nonconvex, we use the iteratively reweighted method for resolving the nonconvex subproblem.
4.1. The ALM Method
Through the introduction of an auxiliary variable , we transform formula (13) to a constrained form:
In the following, to compute the minimizers of formula (15), we adopt the ALM algorithm because of its fast convergence applicability for many problems; the “allowance” variable is merged it into the above objective function, which leads to solving the following minimization problem:
The parameter of formulation (16) does not need to be set large enough; it demonstrates that the ALM method is more robust to the choice of parameters. With regard to our nonconvex case, although there is no strict proof to guarantee the convergence of our iteration process, numerical experiments in Section 5 will demonstrate that our method can produce better denoising results than nonconvex TV regularizer-based method and TGV regularizer-based method.
To solve the minimization problem (16), we use the alternating iteration algorithm. Firstly, we fix and to solve , and then we fix to solve and . More precisely, the alternating iteration algorithm is equivalent to iteratively solving the following two subproblems:
Let the first-order derivative of formula (17) be equal to zero, leading to the following equation:
The minimization problem can be solved by using the Newton iteration of (19). For solving the nonconvex minimization problem (18), we propose the iteratively reweighted method.
4.2. The Iteratively Reweighted Method
Originally, the iteratively reweighted algorithm proceeds by iteratively solving problems which approximate the original problem [17, 20, 21], which was proposed to improve the sparsity in regularized compressed sensing problems, but it turns out that this algorithm is also useful for computer vision applications.
For solving the optimization problem (18), we transform it into the following equivalent energy minimization problem:
As the inner problem (20) is a convex minimization problem, it can be solved efficiently by the primal-dual algorithms [22, 23]. In particular, primal-dual form of problem (20) is derived by the duality principles: where
The primal-dual algorithm for problem (21) reads as follows.
Algorithm 1. The primal-dual algorithm for nonconvex TGV denoising is as follows.(1)Initialization
While do(2)(3)(4) (5) end while.
Set the following:
and is negative conjugate of , .
For the (outer) nonconvex problem, let be the sequence generated by the problem (20), where the index refers to the inner iterations for solving the convex problem and to the outer iterations. According to Proposition 1 [18], will be monotonically decreasing and provides a natural stopping criterion for the inner and outer problem.
The inner iterations are stopped when , where is the maximal number of inner iterations; we set . Then, we stop the outer iterations as soon as
where is a threshold indicating the required accuracy. In our experiments, we use .
Integrating the ALM process and the iteratively reweighted iteration, the whole speckle noise removal algorithm can be found in Algorithm 2.
Algorithm 2. The nonconvex TGV regularizer-based speckle noise removal algorithm is as follows.(1)Initialization is as follows: give parameters ; initialize , .(2)Repeat.(3)Apply Newton method to solve for by formula (19).(4)Apply primal-dual method to solve for by Algorithm 1.(5).(6)Until a stopping criterion is satisfied, final restored image .
5. Numerical Results
In this section, we report numerical experiments on speckle reduction to validate the efficiency of the proposed method with 10 test images whose sizes are all 256*256 shown in Figure 1, including six “standard” test images and four some synthetic images used in [13, 17, 24]. Finally, we test three real SAR images used in [13]. During the evaluation, we compare the proposed algorithm with two recent speckle noise removal algorithms which include NRSNR algorithms [17] and TGVSNR algorithms [13].
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
The evaluating indicator is SNR which is measured between noise-free and restored noisy images. If and are the estimated and noise-free images and is the average variance of the noise-free image , the SNR is defined as where MSE is the mean-square error which is given by .
All simulations listed here are implemented in Matlab R2009a on a laptop equipped with 1.60 GHz CPU and 4 G RAM memory.
5.1. Choice of the Regularization Parameters
Some required parameters , and need to be given to begin with our nonconvex TGV regularization algorithm. The parameter is used to approximate to the constraint in mode (15). In all our experiments, we find the best of each image from a search interval [5, 9]. It is worth noting that the impact of parameter is not great; the choices of are proper for most cases. The parameters and are used to make the primal-dual algorithm’s converging. The two parameters are set too small; the algorithm will be too slow; if the two parameters are set too large, the algorithm will not be robust. We find that will yield good results empirically. The parameters in the primal-dual algorithm have a strong impact on the denoising results. The adjustment of these parameters not only is concerned with noise but also depends on images. We have to turn them manually. We search for the best of each image from a 0.5-step interval , and the best of each image from a 0.3-step interval . As two parameters act on the regularization term, we can fix one of them to tune the other one. Generally, the values of these two regularization parameters are set smaller with the decreasing of the number of looks .
5.2. Influence of the Parameters
Parameters of nonconvex TGV regularization term have an important influence on the performance of the proposed algorithm. Experimental results caused by the different parameters are not the same. We test many images and obtain similar experimental result as in Figure 2. The curves in Figure 2 are based on “House” image. Final SNR of restoration images firstly increases and then decreases with the increasing of value shown in Figure 2, so value for all the experiments is set to .
5.3. Results with Known Image
In this section, we measure the despeckling capability on experimental speckled images whose reference originals are known. Table 1 presents the comparison with only two recent speckle noise removal methods, because these two recent algorithms which include NRSNR algorithms [17] and TGVSNR algorithms [13] are superior to many results before obtained including the algorithms of the AA method and the BF method. The best SNR of different algorithms are listed in Table 1. From the data in the table, we can find that for almost all the denoising results the proposed algorithm achieves higher values of SNR (averagely exceeding about 0.49 db over the TGVSNR method and 0.53 db over the NRSNR method). Even for the “worst” result, our algorithm yields comparable SNR compared with the best values obtained from the latter two methods. Consequently, we believe that the proposed method can averagely perform better than the other two methods. In Figures 3 and 4, we display four representative denoising results to manifest denoising effect and the edge-preserving capability of the proposed method. From observation, we can find that the images restored from the proposed method (shown in Figure 3(b)) are better than those restored from the TGVSNR method (shown in Figure 3(c)) and the NRSNR method (shown in Figure 3(d)). As a matter of fact, in Figure 3(d), we observe that severe staircasing artifacts occur in the restored images of the NRSNR method. In some ways, the TGVSNR method generates better denoising results than the NRSNR method; particularly, staircasing artifacts are weakened greatly. However, some slight residual noise is visible in the restored images of the TGVSNR method. Compared with the NRSNR method and the TGVSNR method, the proposed method generates better denoising results; not only is the speckle noise well cleaned but also the details of the images are well preserved. The above arguments can be also verified from the amplification of the despeckled images shown in Figure 4. From the local zoomed results, we can take a closer look at the advantage of our method more directly. The proposed method can remove the staircasing artifacts better which exist in the areas of the NRSNR method compared to TGVSNR method; meanwhile, it maintains better denoising capability and neat edges.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
5.4. Results on Real SAR Images
In this subsection, we test three real SAR images to observe the performance of the proposed method. As the clean image is unknown, the corresponding algorithms will be terminated when the algorithm reaches the stopping criterion. All of the three methods remove the speckle in homogeneous regions and preserve edges well. By observation, the nonconvex TV-based NRSNR method result as the last column in Figure 5 introduces visible staircasing artifacts; we can find that the details of the images restored by our method are clearer than those of the TGVSNR method and the NRSNR method. It further demonstrates that the proposed nonconvex TGV regularizer-based method can better preserve details of images compared to the other two methods.
(a)
(b)
(c)
(d)
6. Conclusions
In this paper, we propose a novel despeckling variation method based on nonconvex TGV regularization, and an efficient solving algorithm originated from the ALM iteration and the iteratively reweighted primal-dual method is used in our method. Numerical experiments demonstrate that the proposed denoising method based on the nonconvex TGV regularizer shows better performance than nonconvex TV regularizer-based method and TGV regularizer-based method in the ways of preserving the details and edges of images.
In addition, our proposed algorithm is not perfect; one drawback is that its running time is longer than the other compared algorithms, so more rapid algorithm is necessary to be hope for in the future.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
The authors would like to thank Dr. Wensen Feng and Yu Han for sharing their paper codes with them. This paper is supported by the NSFC 61271294 and 61362029 and Shaanxi Provincial Department of Education Project 14JK1827.