Second-order TGV model for Poisson noise image restoration.

Restoring Poissonian noise images have drawn a lot of attention in recent years. There are many regularization methods to solve this problem and one of the most famous methods is the total variation model. In this paper, by adding a quadratic regularization on TGV regularization part, a new image restoration model is proposed based on second-order total generalized variation regularization. Then the split Bregman iteration algorithm was used to solve this new model. The experimental results show that the proposed model and algorithm can deal with Poisson image restoration problem well. What's more, the restoration model performance is significantly improved both in visual effect and objective evaluation indexes.

Background

Image restoration is a fundamental and unavoidable problem in digital image processing by reconstructing the blurred or noisy image from the observed image and not losing the details of image (Zhao et al. 2013; Bilal et al. 2015). There are abundant methods to restore these images. For example, in Zhao et al. (2013), a Bayesian minimum mean squared error (MMSE) estimation based on high-order non-local range-Markov random field (NLRMRF) prior is proposed for non-blind image deblurring problem. The constrained optimization framework (Bilal et al. 2015) was presented to solve the image spatial degradation problem. Moreover, an adaptive weighted regularization scheme was also proposed in modified error estimate (MEE) to cater with the uncertainty due to ill-posed nature of the inverse problem in Bilal et al. (2015). In terms of regularization models, one of the most remarkable models was the ROF model introduced by Rudin et al. (1992), which presented total variation regularization for removing the additive Gaussian white noise. The total variation (TV) model assumed the original image u to be defined on and obtained a solution from a minimization problemwhere the first term is the total variation regularization term, the second term is fidelity term and is a positive fidelity parameter. Although the TV model is effective in image processing, specially for Gaussian white noise. However, this model is not effective in restoring Poisson noise images such as astronomical (Starck and Murtagh 2007), biomedical (Dey et al. 2006; Sarder and Nehorai 2006; Willett and Nowak 2003), and photographic imaging (Foi et al. 2005). For this reason, scholars have proposed some new methods to deal with Poisson noise image problems.

Le et al. (2007) proposed a new total variation model to recover image corrupted by Poisson noise, the new total variation model with fidelity term is suitable for Poisson noise. The new model can be written asThe authors used the gradient descent method to obtain its optimum solution. However, this method can not obtain optimal approximation when the image is both high intensity noise and low intensity features. In Sawatzky et al. (2013), an efficient EM-TV algorithm is presented to speed the computation of the optimization problem (2). In addition, alternating split Bregman iterative algorithm (Setzer et al. 2010) is also used to solve the question (2), since this algorithm does not contain iterations and also not produce negative values. In addition, Figueiredo and Bioucas-Dias (2010) proposed an approach based on alternating direction optimization method for deconvolving Poissonian images.

Recently, based on Chavent and Kunisch (1997) and Liu and Huang (2012) proposed another new total bounded variation-based Poissonian images restoration modelAlthough the above model is better than the total variation model and has a competitive superiority, there also exist some shortcomings, for example, sometimes it will cause undesired oil painting artifacts. In order to avoid the staircase effect, many methods have been proposed. A well-known method to eliminate staircase effect is the TGV (Bredies et al. 2010; Bredies and Valkonen 2011; Bredies et al. 2013) regularization. The TGV regularizer can effectively eliminate the staircase effect but there are still some shortcomings, it tends to introduce some blurring on image edges and texture regions as the existence of high-order derivative term. More seriously, some small details will be lost during the denoising. For this reason, we consider combine the TGV and as one regularization term to solve the poisson noise image restoration problem.

The rest of this article is organized as follows. In “Total generalized variation (TGV)” section, we briefly review the total generalized variation (TGV). The proposed model and algorithm are presented in “The proposed Poisson noise recovering model and algorithm” section. In “Experimental results and discussions” section, experimental results are illustrated to show the consistent performance of the proposed method. Finally, conclusions are given in “Conclusions” section.

Total generalized variation (TGV)

Bredies et al. (2010) proposed the concept of total generalized variation (TGV), which is considered to be the generalization of TV. For convenience, some concepts of TGV are given as follows.

Definition 1

(Bredies et al. 2010) Let be a domain, and . Then, the total generalized variation of order k with weight for is defined as the value of the functionwhere denotes the space of symmetric tensors of order k with arguments in , and are fixed positive parameters.

Definition 2

(Bredies et al. 2010) The space of bounded generalized variation is defined asHere is a Banach space independent of the weight vector .

Definition 3

(Bredies et al. 2010) The “dualization” in the definition of the functional can also be informally interpreted in terms of iterated Fenchel duality.Note that the tensor field are in different spaces for varying l. Moreover, the operator denotes the symmetrized gradient operator

In this paper, we use in the proposed model. Thus, the second-order TGV can be written aswhere denotes the space of symmetric matrices. And the first and second divergences are defined asIn addition, according to Bredies and Valkonen (2011), the energy term can be formulated aswhere can be separately expressed as

The proposed Poisson noise recovering model and algorithm

The proposed model for Poisson noise image

We assume that is the original image, is the observed image, is a linear blurring operator related with the spread point function (PSF). Then the degradation model can be described aswhere P denotes the Poisson distribution function.

Based on Le et al. (2007), we define the Bayes Law as follows.According to (14), for each , we haveNext, we assume that the prior distribution P(u) is TGV and , which can be written aswhere is the regularization parameter.

Thus, we obtain a model for restoring the Poissonian noise image as follows.where is the indicator function of set SBy reformulating TGV as a minimization in the discrete setting, the proposed model can be written as

The split Bregman algorithm for Poisson noise removal

The split Bregman algorithm (Goldstein and Osher 2009; Wang et al. 2008) has been widely used in image processing, which is easy to be realized and has fast convergence (Cai et al. 2009; Jia et al. 2009). Therefore, we use split Bregman algorithm to solve our minimization problem (19). Firstly, by introducing new auxiliary variables w, x, y and z, the problem (19) can be reformulated as the following constrained optimization problemFor the above constrained problem (20), we transform it into the corresponding unconstrained problemwhere are positive penalty parameters. Thus the split Bregma iterative algorithm for solving the question (20) can be described aswhere the updates of the multipliers is described as followsSince the updates of are merely simple calculations, then the minimization question (22) can be divided into the following several subproblems: Given initial value and , the split Bregman algorithm can be written asFor w-subproblem, note that it is separable with respect to each component. It is easy to solve and the solution of w may be written asAs for solving x, y-subproblem, we can directly obtain the solutions by using shrinkage operator:

The x-subproblem can be solved byThe solution of the y-subproblem is similarly obtainedThe (u, p)-subproblem is a saddle-point problem, which can be divided into the following two subproblems:

For u, we havewhich can be solved by considering the following normal equationFinally, u is solved by

For the sub-problem p, it can be written as the following minimization problemwhere is a vector, is a matrix.

For , it can be solved by considering the following linear systemTherefore,Similarly, we can obtain the solution of as

Experimental results and discussions

In this section, we illustrate some numerical results of the proposed model for the Poisson noise removal problem. We compare our method with the one proposed in Figueiredo and Bioucas-Dias (2010) (PIDAL) and the other proposed in Liu and Huang (2012) (PID-Split). In order to prove the superiority of the proposed model, we compare our model with TGV regularization model.To show the effectivity of the proposed model, we choose four pictures possed abundant detail information.

We terminate the iterations for these methods by the following stopping criterionThe quality of the restoration results is compared quantitatively by using the Signal-to-Noise Ratio (SNR), the Peak Signal-to-Noise (PSNR), the relative error (RelErr) and the Structural SIMilarity index (SSIM). They are defined as followswhere u and are the ideal image and the restored image, respectively.where and are averages of u and , respectively. and are the variance of u and , respectively. is the covariance of u and . The positive constants and can be thought of as stabilizing constants for near-zero denominator values. Generally speaking, the more bigger value of SNR, PSNR or the smaller value of RelErr is, the better quality of the reconstructed image is.

The Poissonian images used for our experiments are generated as follows: the original images are convoluted with a blur kernel and additionally contaminated by Poisson noise, here we use the function in MATLAB’s Statistics Toolbox after blurring the true images with the given point spread functions to generate the blurred and noise images.

The selection of the regularization parameters highly affects the image restoration results, and related to make the fair comparison with different denoising models. The penalty parameters which relies on unknown noise level highly influences the speed of the algorithms. In experiments, we set in the PIDAL algorithm. In the PID-Split algorithms, we choose . In the TGV model, we set . The penalty parameter in the proposed method is empirically set . Thus, we may have a good restoration results.

In the first experiment, we used the image “Woman” () in Fig. 1a. We perform the blurring operation psfGauss(5, 2) proposed in Nagy et al. (2004) on the original image and add the Poisson noise to the blurred data to generate the degraded image in Fig. 1b. The parameter of this test, we set in PIDAL algorithm, for PID-Split algorithms,due to the TGV model we set , set for the proposed model. The pictures of Fig. 1c–f are the restoration images, which represent the difference between the three methods. From these pictures, we can see the proposed model have more advantages. In order to more effectively reflect the experiment result, Fig. 1g–j present the residual images refer to the difference of the original image and the restoration image. From these pictures, we can see that the proposed model can preserve more details than other methods. In the Table 1, the SNR, PSNR, RelErr and SSIM values of the restored images by the proposed model are better than other methods.

Fig. 1

The woman picture is compared with other method. a Original image; b degraded image; c the PIDAL model result; d the PID-Split method result; e the TGV method result; f the proposed model result. g–j The residual images

Table 1

Summarized all of the experiment restoration results

Method	SNR	PSNR	RelErr	SSIM
Test1
PIDAL (Figueiredo and Bioucas-Dias 2010)	25.33	30.30	0.054	0.904
PID-Split (Liu and Huang 2012)	25.18	30.16	0.055	0.917
TGV	25.31	30.32	0.053	0.908
Proposed	25.48	30.45	0.053	0.922
Test2
PIDAL (Figueiredo and Bioucas-Dias 2010)	21.17	25.03	0.087	0.828
PID-Split (Liu and Huang 2012)	21.38	25.36	0.084	0.837
TGV	21.24	25.10	0.086	0.835
Proposed	21.74	25.60	0.080	0.843
Test3
PIDAL (Figueiredo and Bioucas-Dias 2010)	21.50	27.87	0.0841	0.811
PID-Split (Liu and Huang 2012)	21.49	27.86	0.0357	0.811
TGV	20.64	27.01	0.092	0.801
Proposed	21.01	27.23	0.0803	0.820

In second experiment, we use the image “House” with size () in Fig. 2a, which contains a lot of edge details. We perform the blurring operation with radius 5 () on the original image and add the Poisson noise to the blurred data to generate the degraded image in Fig. 2b. As for parameter selection, we choose for the PIDAL algorithm, for the PID-Split algorithm, set in the TGV model. The proposed algorithm we set . Form Fig. 2c–f, we can see that the proposed model compared to the PIDAL method and PID-Split algorithms have better restoration results. In Fig. 2g–j, we have enlarged some details of the images, which can be clearly see the advantages of the proposed model for the recovery of edge details. The SNR, RelRrr and SSIM values in Table 1 showed that the proposed model have a better restoration result.

Fig. 2

The house picture is compared with other method. a Original image; b degraded image; c the PIDAL model result; d the PID-Split method result; e the TGV method result; f the proposed model result; g–j The result of partial enlarged pictures

In order to further verify the validity of the new model, in third experiment, we use the image which contains a lot of edge details. We perform the blurring operation by a line motion blur. The point spread function for the linear motion blur is returns a filter to approximate, once convolved with an image, the linear motion of a camera by r pixels, with an angle of degrees in a counter-clockwise direction. In this example, and , then add the Poisson noise to the blurred data to generate the degraded image in Fig. 3b. The parameters choose as the same as those in second experiment and also may be adjusted.

Fig. 3

The third experiment. a Original image; b degraded image; c the PIDAL model result; d the PID-Split method result; e the TGV method result; f the proposed model result; g–j The residual images

Finally, let us choose a human brain MR image of size as the test image. We zoom in a marked close-up region which is abundant in texture-like features to better visual comparison. We can clearly see that the produced textures by our proposed method are better quality than the other methods from Fig. 4.

Fig. 4

The fourth experiment. a Ground truth; b ground truth enlarge; c PID-Split method result; d PID-Split method enlarge; e PIDAL-Method result; f PIDAL-method enlarge; g Proposed method result; h proposed method enlarge

Conclusions

In this paper, we investigate the second-order total generalized variation with a quadratic regularization to deal with the Poissonian images restoration problem. The proposed model is solved efficiently by split Bregman iterative algorithm in this way the calculation speed is fast. Numerical results show that our proposed method is particularly advantageous for restoration the Poisson images in terms of SNR, SSIM and RelErr quality compared to other methods. In the model, the parameters selection is a difficult problem which needs further study.