Image Restoration by Iterative Denoising and Backward Projections
Abstract
Inverse problems appear in many applications such as image deblurring and inpainting. The common approach to address them is to design a specific algorithm for each problem. The PlugandPlay (P&P) framework, which has been recently introduced, allows solving general inverse problems by leveraging the impressive capabilities of existing denoising algorithms. While this fresh strategy has found many applications, a burdensome parameter tuning is often required in order to obtain highquality results. In this work, we propose an alternative method for solving inverse problems using denoising algorithms, that requires less parameter tuning. We demonstrate that it is competitive with taskspecific techniques and the P&P approach for image inpainting and deblurring.
I Introduction
We consider the reconstruction of an image from its degraded version, which may be noisy, blurred, downsampled, or all together. This general problem has many important applications, such as medical imaging, surveillance, entertainment, and more. Traditionally, the design of taskspecific algorithms has been the ruling approach. Many works specifically considered the denoising problem [1, 2, 3], the deblurring problem [4, 5, 6], the inpainting problem [7, 8], etc.
Recently, a new approach attracts much interest. This approach suggests leveraging the impressive capabilities of existing denoising algorithms for solving other tasks that can be formulated as an inverse problem. The pioneering algorithm that introduced this concept is the PlugandPlay (P&P) method [9], which presents an elegant way to decouple the measurement model and the image prior, such that the latter is handled solely by a denoising operation. Thus, it is not required to explicitly specify the prior, since it is implicitly defined through the choice of the denoiser.
The P&P method has already found many applications, e.g. bright field electron tomography [10], Poisson denoising [11], and postprocessing of compressed images [12]. It also inspired new related techniques [13, 14, 15]. However, it has been noticed that the P&P often requires a burdensome parameter tuning in order to obtain high quality results [14, 16]. Moreover, since it is an iterative method, sometimes a large number of iterations is required.
In this work, we propose a simple iterative method for solving linear inverse problems using denoising algorithms, which provides an alternative to P&P. Our strategy has less parameters that require tuning (e.g. no tuning is required for the noisy inpainting problem), often requires less iterations, and its recovery performance is competitive with taskspecific algorithms and with the P&P approach. We demonstrate the advantages of the new technique on inpainting and deblurring problems.
The paper is organized as follows. In Section II we present the problem formulation and the P&P approach. The proposed algorithm is presented in Section III. Section IV includes mathematical analysis of the algorithm and provides a practical way to tune its parameter. In Section V the usage of the method is demonstrated and examined for inpainting and deblurring problems. Section VI concludes the paper.
Ii Background
Iia Problem formulation
The problem of image restoration can be generally formulated by
(1) 
where represents the unknown original image, represents the observations, is an degradation matrix and is a vector of independent and identically distributed Gaussian random variables with zero mean and standard deviation of . The model in (1) can represent different image restoration problems; for example: image denoising when is the identity matrix , image inpainting when is a selection of rows of , and image deblurring when is a blurring operator.
In all of these cases, a prior image model is required in order to successfully estimate from the observations . Specifically, note that is illconditioned in the case of image deblurring, thus, in practice it can be approximated by a rankdeficient matrix, or alternatively by a full rank matrix (). Therefore, for a unified formulation of inpainting and deblurring problems, which are the test cases of this paper, we assume .
Almost any approach for recovering involves formulating a cost function, composed of fidelity and penalty terms, which is minimized by the desired solution. The fidelity term ensures that the solution agrees with the measurements, and is often derived from the negative loglikelihood function. The penalty term regularizes the optimization problem through the prior image model . Hence, the typical cost function is
(2) 
where stands for the Euclidean norm.
IiB Plug and Play approach
Instead of devising a separate algorithm to solve for each type of matrix , a general recovery strategy has been proposed in [9], denoted as the PlugandPlay (P&P). For completeness, we briefly describe this technique.
Using variable splitting, the P&P method restates the minimization problem as
(3) 
where is the fidelity term in (2), and is a positive parameter that adds flexibility to the cost function. This problem can be solved using ADMM [17] by constructing an augmented Lagrangian, which is given by
(4) 
where is the dual variable, is the scaled dual variable, and is the ADMM penalty parameter. The ADMM algorithm consists of iterating until convergence over the following three steps
(5) 
By plugging (IIB) in (IIB) we have
(6) 
Note that the first step in (IIB) is just solving a least squares (LS) problem and the third step is a simple update. The second step is more interesting. It describes obtaining using a white Gaussian denoiser with noise variance of , applied on the image . This can be written compactly as , where is a denoising operator. Since general denoising algorithms can be used to implement the operator , the P&P method does not require knowing or explicitly specifying the prior function . Instead, is implicitly defined through the choice of . The obtained P&P algorithm is presented in Algorithm 1.
The convergence of the P&P method is proved for Gaussian denoisers that satisfy certain conditions [10]. However, well known denoisers such as BM3D [1], KSVD [2], and standard NLM [3], lead to good results despite violating these conditions.
The P&P method is not free of drawbacks. Its main difficulties are the large number of iterations, which is often required by the P&P to converge to a good solution, and the setting of the design parameters and , which is not always clear and strongly affects the performance.
Iii The Proposed Algorithm
In this work we take another strategy for solving inverse problems using denoising algorithms. We start with formulating the cost function (2) in somewhat strange but equivalent way
(7) 
where
(8)  
(9) 
Note that is the pseudoinverse of the full row rank matrix , and is not a real norm, since is not a positive definite matrix in our case. Moreover, as mentioned above, since the null space of is not empty, the prior is essential in order to obtain a meaningful solution.
The optimization problem can be equivalently written as
(10) 
Note that due to the degenerate constraint, the solution for is trivial .
Now, we make two major modifications to the above optimization problem. The basic idea is to loose the variable in a restricted manner, that can facilitate the estimation of . First, we give some degrees of freedom to by using the constraint instead of . Next, we turn to prevent large components of in the null space of that may strongly disagree with the prior . We do it by replacing the multiplication by in the fidelity term, which implies projection onto a subspace, with multiplication by that implies a full dimensional space, where is a design parameter.
This leads to the following optimization problem
(11) 
Note that introduces a tradeoff. On the one hand, exaggerated value of should be avoided, as it may overreduce the effect of the fidelity term. On the other hand, too small value of may overpenalize unless it is very close to the affine subspace . This limits the effective feasible set of in problem (11), such that it may not include potential solutions of the original problem (10). Therefore, we suggest setting the value of as
(12) 
where denotes the feasible set of problem (11). Note that the feasibility of is dictated by and the feasibility of is dictated by the constraint in (11). The problem of obtaining such value for (or an approximation) is discussed in Section IVA, where a relaxed version of the condition in (III) is presented.
Assuming that solves (III), the property that for feasible and , together with the fact that is one of the solutions of the underdetermined system , prevents increasing the penalty on potential solutions of the original optimization problem (10). Therefore, roughly speaking, we do not lose solutions when we solve (11) instead of (10).
We solve (11) using alternating minimization. Iteratively, is estimated by solving
(13) 
and is estimated by solving
(14) 
which describes a projection of onto the affine subspace , and has a closedform solution
(15) 
Similarly to the P&P technique, (13) describes obtaining using a white Gaussian denoiser with noise variance of , applied on the image , and can be written compactly as , where is a denoising operator. Moreover, as in the case of the P&P, the proposed method does not require knowing or explicitly specifying the prior function . Instead, is implicitly defined through the choice of .
The variable is expected to be closer to the true signal than the raw observations . Thus, our algorithm alternates between estimating the signal and using this estimation in order to obtain improved measurements (that also comply with the original observations ). The proposed algorithm, which we call Iterative Denoising and Backward Projections (IDBP), is presented in Algorithm 2.
Iv Mathematical Analysis of the Algorithm
Iva Setting the value of the parameter
Setting the value of that solves (III) is required for simple theoretical justification of our method. However, it is not clear how to obtain such in general. Therefore, in order to relax the condition in (III), that should be satisfied by all and in , we can focus only on the sequences and generated by the proposed alternating minimization process. Then, we can use the following proposition.
Proposition 1.
Set . If there exist an iteration of IDBP that violates the following condition
(16) 
then also violates the condition in (III).
Proof.
Assume that and generated by IDBP at some iteration violate (16), then they also violate the equivalent condition
(17) 
Note that (17) is obtained simply by plugging (15) into in (III). Therefore, and also violate the inequality in (III). Finally, it is easy to see that and are feasible points of (11), since is a feasible point of and satisfies . Therefore, the condition in (III) does not hold for all feasible and , which means that violates it. ∎
Note that (16) can be easily evaluated for each iteration. Thus, violation of (III) can be spotted (by violation of (16)) and used for stopping the process, increasing and running the algorithm again. Of course, the opposite direction does not hold. Even when (16) is satisfied for all iterations, it does not guarantee satisfying (III). However, the relaxed condition (16) provides an easy way to set with an approximation to the solution of (III), which gives very good results in our experiments.
In the special case of the inpainting problem, (16) becomes ridiculously simple. Since is a selection of rows of , it follows that , which is an matrix that merely pads with zeros the vector on which it is applied. Therefore, , implying that satisfies (16) in this case. Obviously, if , a small positive is required in order to prevent the algorithm from getting stuck (because in this case ).
Condition (16) is more complex when considering the deblurring problem. In this case is an illconditioned matrix. Therefore must be approximated, either by approximating by a full rank matrix before computing (8), or by regularized inversion techniques for , e.g. standard Tikhonov regularization. A research on how to compute in this case is ongoing. We empirically observed that using a fixed value for (for all noise levels and blur kernels) exhibits good performance. However, we had to add another parameter that controls the amount of regularization in the approximation of , that we slightly change between scenarios (i.e. when or change). This issue is discussed in Section VB. An interesting observation is that the pairs of which give the best results indeed satisfy condition (16). On the other hand, the pairs of that give bad results often violate this condition (recall that the condition should be met during all iterations). An example of this behavior is given in the end of Section VB.
IvB Analysis of the sequence
The IDBP algorithm creates the sequence that can be interpreted as a sequence of updated measurements. It is desired that is improved with each iteration, i.e. that , obtained from , estimates better than , which is obtained from .
Assuming that the result of the denoiser, denoted by , is perfect, i.e. , we get from (15)
(18) 
The last equality describes a model that has only noise (possibly colored), and is much easier to deal with than the original model (1). Therefore, can be considered as the optimal improved measurements that our algorithm can achieve. Since we wish to make no specific assumptions on the denoising scheme , improvement of will be measured by the Euclidean distance to .
Let us define the orthogonal projection onto the row space of , and its orthogonal complement . The updated measurements are always consistent with on , and do not depend on , as can be seen from
(19) 
Thus, the following theorem ensures that iteration improves the results, provided that is closer to than on the null space of , i.e.
(20) 
Theorem 2.
Assuming that (20) holds at the th iteration of IDBP, then we have
(21) 
Proof.
A denoiser that makes use of a good prior (and suitable ) is expected to satisfy (20), at least in early iterations. For example, in the inpainting problem is associated with the missing pixels, and in the deblurring problem is associated with the data that suffer the greatest loss by the blur kernel. Therefore, in both cases is expected to be closer to than . Note that if (20) holds for all iterations, then Theorem 2 ensures monotonically improvement and convergence of . However, it still does not guarantee that is the limit of the sequence.
V Experiments
We demonstrate the usage of IDBP for two test scenarios: the inpainting and the deblurring problems. We compare the IDBP performance to P&P and another algorithm that has been specially tailored for each problem [6], [18]. In all experiments we use BM3D [1] as the denoising algorithm for IDBP and P&P. We use the following four test images in all experiments: cameraman, house, peppers and Lena. Their intensity range is 0255.
Va Image inpainting
In the image inpainting problem, is a selection of rows of and , which simplifies both P&P and IDBP. In P&P, the first step can be solved for each pixel individually. In IDBP, is obtained merely by taking the observed pixels from and the missing pixels from . For both methods we use the result of a simple median scheme as their initialization (for in P&P and for in IDBP). It is also possible to alternatively use for initialization, but then many more iterations are required. Note that the computational cost of each iteration of P&P and IDBP is of the same scale, dominated by the complexity of the denoising operation.
The first experiment demonstrates the performance of IDBP, P&P and inpainting based on Image Processing using Patch Ordering (IPPO) approach [18], for the noiseless case () with 80% missing pixels, selected at random. The parameters of IPPO are set exactly as in [18], where the same scenario is examined. The parameters of P&P are optimized for best reconstruction quality. We use , and 150 iterations. Also, for P&P we assume that the noise standard deviation is 0.001, i.e. nonzero, in order to compute .
Considering IDBP, in Section IVA, it is suggested that . However, since in this case , a small positive , e.g. , is required. Indeed, this setting gives good performance, but also requires ten times more iterations than P&P. Therefore, we use an alternative approach. We set , which allows us to use only 150 iterations (same as P&P), but take the last as the final estimate, which is equivalent to performing the last denoising with the recommended . Figure 1 shows the results of both IDBP implementations for the house image. It approves that the alternative implementation performs well and requires significantly less iterations (note that the xaxis has a logarithmic scale). Therefore, for the comparison of the different inpainting methods in this experiment, we use IDBP with .
The results of the three algorithms are given in Table I. IDBP is usually better than IPPO, but slightly inferior to P&P. This is the cost of enhancing IDBP by setting to a value which is significantly larger than zero. However, this observation also hints that IDBP may shine for noisy measurements, where can be used without increasing the number of iterations. We also remark that IPPO gives the best result for peppers because in this image P&P and IDBP require more than the fixed 150 iterations.
camera.  house  peppers  Lena  
IPPO  24.78  32.64  27.98  31.84 
P&P  24.83  34.72  26.88  32.41 
IDBP  24.86  33.78  26.86  32.13 
The second experiment demonstrates the performance of IDBP and P&P with 80% missing pixels, as before, but this time . Noisy inpainting has not been implemented yet by IPPO [18]. The parameters of P&P that give us the best results are , and 150 iterations. Using the same parameter values as before deteriorates the performance significantly. Contrary to P&P, in this experiment tuning the parameters of IDBP can be avoided. We follow Section IVA and set . Moreover, IDBP now requires only 75 iterations, half the number of P&P. The results are given in Table II. P&P is slightly inferior to IDBP, despite having twice the number of iterations and a burdensome parameter tuning. The results for house are also presented in Figure 2, where it can be seen that P&P reconstruction suffers from more artifacts (e.g. ringing artifacts near the right window).
camera.  house  peppers  Lena  

P&P  24.55  31.53  26.16  30.10 
IDBP  24.68  31.62  26.13  30.14 
We repeat the last experiment with slightly increased noise level of , but still use the same parameter tuning for P&P, which is optimized for (i.e. , and the fixed ). This situation is often encountered in practice, when calibrating a system for all possible scenarios is impossible. The results are given in Table III. The IDBP clearly outperforms P&P in this case. This experiment clearly shows the main advantage of our algorithm over P&P as it is less sensitive to parameter tuning.
camera.  house  peppers  Lena  

P&P  24.43  30.78  25.80  29.47 
IDBP  24.51  31.14  25.92  29.69 
VB Image deblurring
In the image deblurring problem, for a circular shiftinvariant blur operator whose kernel is , both P&P and IDBP can be efficiently implemented using Fast Fourier Transform (FFT). In P&P, can be computed by
(24) 
where denotes the FFT operator, denotes the inverse FFT operator.
Recall that is an illconditioned matrix. Therefore, In IDBP we replace with a regularized inversion of , using standard Tikhonov regularization, which is given in the Fourier domain by
(25) 
where is a parameter that controls the amount of regularization in the approximation of . Then, in IDBP can be computed by
(26) 
While a research on how to compute in this case is ongoing, we empirically observed that using a fixed value exhibits good performance with only slightly changing between the examined scenarios (i.e. different configurations of and ). Lastly, we remark that we use trivial initialization in both methods, i.e. in P&P and in IDBP. Similarly to inpainting, the computational cost of each iteration of P&P and IDBP is on the same scale, dominated by the complexity of the denoising operation.
We consider four deblurring scenarios used as benchmarks in many publications (e.g. [5, 6]). The blur kernel and noise level of each scenario are summarized in Table IV. The kernels are normalized such that . Table V shows the results of IDBP, P&P and the dedicated algorithm BM3DDEB [6]. For each scenario it shows the input PSNR (i.e. PSNR of ) and the BSNR (blurred signaltonoiseratio, defined as ) for each image, as well as the ISNR (improvement signaltonoiseratio) for each method and image, which is the difference between the PSNR of the reconstruction and the input PSNR. Note that in Scenario 3, is set slightly different for each image, ensuring that the BSNR is 40 dB.
Scenario  

1  2  
2  8  
3  uniform  0.3 
4  49 
The parameters of BM3DDEB are set exactly as in [6], where the same scenarios are examined. The parameters of P&P are optimized for each scenario. We have 0.85, 0.85, 0.9, 0.8 and 2, 3, 3, 1, for scenarios 14, respectively, and use 50 iterations. IDBP is easier to tune, as we set 7e3, 4e3, 8e3, 2e3 for scenarios 14, respectively, and fixed . Also, we use only 20 iterations for IDBP. From Table V, it is clear that IDBP and P&P are highly competitive (arguably, IDBP performs slightly better), and both are much better than BM3DDEB, which is tailored for the deblurring problem. Figure 3 displays the results for Lena in Scenario 2. It can be seen that IDBP reconstruction is the sharpest (e.g. notice the lightning in the right eye).
As mentioned in Section IVA, we observed that the pairs of which give the best results indeed satisfy condition (16), while the pairs of that give bad results often violate this condition. This behavior is demonstrated for house image in Scenario 1. Figure 3(a) shows the PSNR vs. iteration for several pairs of . The lefthand side of (16) divided by the righthand side is presented in Figure 3(b). If this division is less than 1, even for a single iteration, it means that the original condition in (III) is violated by the associated . Recall that even when the division is higher than 1 for all iterations, it does not guarantee satisfying (III). Therefore, a small margin should be kept.
Scenario 1  camera.  house  peppers  Lena 
BSNR  31.87  29.16  29.99  29.89 
input PSNR  22.23  25.61  22.60  27.25 
BM3DDEB  8.20  9.34  9.58  7.96 
P&P  8.03  9.74  10.02  8.02 
IDBP  8.48  9.82  9.86  8.01 
Scenario 2  camera.  house  peppers  Lena 
BSNR  25.85  23.14  23.97  23.87 
input PSNR  22.16  25.46  22.53  27.04 
BM3DDEB  6.46  8.22  7.88  6.55 
P&P  6.70  8.02  7.48  6.35 
IDBP  6.53  8.19  7.67  6.64 
Scenario 3  camera.  house  peppers  Lena 
BSNR  40.00  40.00  40.00  40.00 
input PSNR  20.77  24.11  21.33  25.84 
BM3DDEB  8.37  10.94  10.05  7.99 
P&P  9.49  13.17  11.70  9.04 
IDBP  9.55  12.90  11.88  9.00 
Scenario 4  camera.  house  peppers  Lena 
BSNR  18.53  15.99  17.01  16.47 
input PSNR  24.63  28.06  24.77  28.81 
BM3DDEB  3.35  5.20  3.31  4.80 
P&P  3.31  5.43  4.95  4.84 
IDBP  3.59  5.66  4.12  5.05 





Vi Conclusion
We presented the Iterative Denoising and Backward Projections (IDBP) method for solving linear inverse problems using denoising algorithms. This method, in its general form, has only a single parameter that should be set according to a given condition. We presented mathematical analysis of this strategy and provided a practical way to tune its parameter. Therefore, it can be argued that our method has less parameters that require tuning than the P&P method, especially for the noisy inpainting problem, where the single parameter of the IDBP can be just set to zero. Experiments demonstrated that IDBP is competitive with taskspecific algorithms and with the P&P approach for inpainting and deblurring problems.
References
 [1] K. Dabov, A. Foi, V. Katkovnik, and K. Egiazarian, “Image denoising by sparse 3D transformdomain collaborative filtering,” IEEE Transactions on image processing, vol. 16, no. 8, pp. 2080–2095, 2007.
 [2] M. Elad and M. Aharon, “Image denoising via sparse and redundant representations over learned dictionaries,” IEEE Transactions on Image processing, vol. 15, no. 12, pp. 3736–3745, 2006.
 [3] A. Buades, B. Coll, and J.M. Morel, “A review of image denoising algorithms, with a new one,” Multiscale Modeling & Simulation, vol. 4, no. 2, pp. 490–530, 2005.
 [4] M. Delbracio and G. Sapiro, “Burst deblurring: Removing camera shake through fourier burst accumulation,” in Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pp. 2385–2393, 2015.
 [5] J. A. GuerreroColón, L. Mancera, and J. Portilla, “Image restoration using spacevariant Gaussian scale mixtures in overcomplete pyramids,” IEEE Transactions on Image Processing, vol. 17, no. 1, pp. 27–41, 2008.
 [6] K. Dabov, A. Foi, V. Katkovnik, and K. O. Egiazarian, “Image restoration by sparse 3D transformdomain collaborative filtering,” in SPIE Electronic Imaging ’08, vol. 6812, (San Jose, California, USA), Jan. 2008.
 [7] M. Bertalmio, G. Sapiro, V. Caselles, and C. Ballester, “Image inpainting,” in Proceedings of the 27th annual conference on Computer graphics and interactive techniques, pp. 417–424, ACM Press/AddisonWesley Publishing Co., 2000.
 [8] A. Criminisi, P. Pérez, and K. Toyama, “Region filling and object removal by exemplarbased image inpainting,” IEEE Transactions on image processing, vol. 13, no. 9, pp. 1200–1212, 2004.
 [9] S. V. Venkatakrishnan, C. A. Bouman, and B. Wohlberg, “Plugandplay priors for model based reconstruction,” in Global Conference on Signal and Information Processing (GlobalSIP), 2013 IEEE, pp. 945–948, IEEE, 2013.
 [10] S. Sreehari, S. V. Venkatakrishnan, B. Wohlberg, G. T. Buzzard, L. F. Drummy, J. P. Simmons, and C. A. Bouman, “Plugandplay priors for bright field electron tomography and sparse interpolation,” IEEE Transactions on Computational Imaging, vol. 2, no. 4, pp. 408–423, 2016.
 [11] A. Rond, R. Giryes, and M. Elad, “Poisson inverse problems by the plugandplay scheme,” Journal of Visual Communication and Image Representation, vol. 41, pp. 96–108, 2016.
 [12] Y. Dar, A. M. Bruckstein, M. Elad, and R. Giryes, “Postprocessing of compressed images via sequential denoising,” IEEE Transactions on Image Processing, vol. 25, no. 7, pp. 3044–3058, 2016.
 [13] T. Meinhardt, M. Möller, C. Hazirbas, and D. Cremers, “Learning proximal operators: Using denoising networks for regularizing inverse imaging problems,” ICCV, 2017.
 [14] Y. Romano, M. Elad, and P. Milanfar, “The little engine that could: Regularization by denoising (red),” arXiv preprint arXiv:1611.02862, 2016.
 [15] A. M. Teodoro, J. M. BioucasDias, and M. A. Figueiredo, “Image restoration and reconstruction using variable splitting and classadapted image priors,” in Image Processing (ICIP), 2016 IEEE International Conference on, pp. 3518–3522, IEEE, 2016.
 [16] S. H. Chan, X. Wang, and O. A. Elgendy, “Plugandplay ADMM for image restoration: Fixedpoint convergence and applications,” IEEE Transactions on Computational Imaging, vol. 3, no. 1, pp. 84–98, 2017.
 [17] S. Boyd, N. Parikh, E. Chu, B. Peleato, and J. Eckstein, “Distributed optimization and statistical learning via the alternating direction method of multipliers,” Foundations and Trends® in Machine Learning, vol. 3, no. 1, pp. 1–122, 2011.
 [18] I. Ram, M. Elad, and I. Cohen, “Image processing using smooth ordering of its patches,” IEEE transactions on image processing, vol. 22, no. 7, pp. 2764–2774, 2013.