External Prior Guided Internal Prior Learning for Real-World Noisy Image Denoising

I Introduction

Image denoising is a crucial and indispensable step to improve image quality in digital imaging systems. In particular, with the decrease of size of CMOS/CCD sensors, image is more easily to be corrupted by noise and hence denoising is becoming increasingly important for high resolution imaging. The problem of image denoising has been extensively studied in literature and numerous image denoising methods have been proposed in the past decades. Most of existing denoising methods focus on the scenario of additive white Gaussian noise (AWGN) , where the observed noisy image $\mathbf{y}$ is modeled as the addition of clean image $\mathbf{x}$ and AWGN $\mathbf{n}$, i.e., $\mathbf{y}=\mathbf{x}+\mathbf{n}$. There are also methods proposed for removing Poisson noise , mixed Poisson and Gaussian noise , mixed Gaussian and impulse noise , and realistic noise in real photography .

Instead of using predefined image priors, methods have also been proposed to learn priors from natural images for denoising. The generative image prior learning methods usually learn prior models from a set of external clean images and apply the learned prior models to the given noisy image , or learn priors from the given noisy image to perform denoising . Recently, the discriminative image prior learning methods , which learn denoising models from pairs of clean and noisy images, have been becoming popular. The representative methods include the neural network based methods , random fields based methods , and reaction diffusion based methods .

Most of the above mentioned methods focus on AWGN removal, however, the assumption of AWGN is too ideal to be true for real-world noisy images, where the noise is much more complex and varies with different scenes, cameras and camera settings (ISO, shutter speed, and aperture, etc.) . As a result, many denoising methods in literature, including those learning based methods, become less effective when applied to real-world noisy images. Fig. 1 shows an example, where we apply some representative and state-of-the-art denoising methods, including CBM3D , WNNM , DnCNN , CSF , and TNRD to a real-world noisy image (captured by a Nikon D800 camera with ISO is 3200) provided in . One can see that these methods either remain much the noise or over-smooth the image details.

There have been a few methods and software toolboxes developed for real-world noisy image denoising. Almost all of these methods follow a two-stage framework: first estimate the parameters of the noise model (usually assumed to be Gaussian or mixture of Gaussians (MoG)), and then perform denoising with the estimated noise model. However, the noise in real-world noisy images is very complex and is hard to be modeled by explicit distributions such as Gaussian and MoG. According to , the noise corrupted in the in-camera imaging process is signal dependent and comes from five main sources: photon shot, fixed pattern, dark current, readout, and quantization noise. The existing methods mentioned above may not perform well on real-world noisy image denoising tasks. Fig. 1 also shows the denoising results of two real-world noisy image denoising methods, Noise Clinic and Neat Image . One can see that these two methods still generate much noise caused artifacts.

This work aims to develop a new paradigm for real-world noisy image denoising. Different from existing real-world noisy image denoising methods which focus on noise modeling, we focus on image prior learning. We argue that with a strong and adaptive prior learning scheme, robust denoising performance on real-world noisy images can still be obtained. To achieve this goal, we propose to first learn image priors from external clean images, and then employ the learned external priors to guide the learning of internal priors from the given noisy image. The flowchart of the proposed method is illustrated in Fig. 2. We first extract millions of patch groups (PGs) from a set of high quality natural images, with which a Gaussian Mixture Model (GMM) is learned as the external image prior. The learned GMM prior model is used to assign each PG extracted from the given noisy image into its most suitable cluster via maximum a-posterior, and then an external-internal hybrid orthogonal dictionary is learned as the final prior for each cluster, with which the denoising can be readily performed by weighted sparse coding with closed form solution. The external priors learned from clean images preserve fine-scale image structural information, which is hard to be reproduced from noisy images. Therefore, external dictionary can serve as a good supplement to the internal dictionary. Our proposed denoising method is simple and efficient, yet our extensive experiments on real-world noisy images demonstrate its better denoising performance than the current state-of-the-arts.

II Related Work

Learning natural image priors plays a key role in image denoising . There are mainly four categories of prior learning based methods. 1) External prior learning methods learn priors (e.g., dictionaries) from a set of external clean images, and the learned priors are used to recover the latent clean image from the given noisy image. 2) Internal prior learning methods directly learn priors from a given noisy image, and image denoising is often done simultaneously with the prior learning process. 3) Discriminative prior learning methods learn discriminative models or mapping functions from clean and noisy image pairs, and the learned models or mapping functions are applied to a noisy image for denoising. 4) Hybrid methods combine the external and internal priors to denoise the given input image.

It has been shown that the external priors learned from natural clean images are effective and efficient for universal image denoising problems, whereas they are not adaptive to the given noisy image and some fine-scale image structures may not be well recovered. By contrast, the internal priors learned from the given noisy image are adaptive to image content, but the learned priors can be much affected by noise and the learning processing is usually slow . Besides, most of the internal prior learning methods assume additive white Gaussian noise (AWGN), making the learned priors less robust for real-world noisy images. In this paper, we use external priors to guide the internal prior learning. Our method is not only much faster than the traditional internal learning methods, but also very robust to denoise real-world noisy images.

In , the authors employed external clean patches to denoise noisy patches with high individual Signal-to-Noise-Ratio (PatchSNR), and employed internal noisy patches to denoise noisy patches with low PatchSNR. This is essentially different from our work which employs the external patch group based prior to guide the clustering and dictionary learning of the internal noisy patch groups. In , the external priors are only used to guide the internal patch clustering for image denoising, while in our work, the learned external priors are employed to guide not only the internal clustering, but also the internal dictionary learning. Besides, the method of follows a patch based framework for AWGN removal, while in our work we employ a patch group based framework for real-world noisy image denoising. In addition, some technical details are also different. For example, method in utilizes low-rank minimization for denoising, while we use dictionary learning and sparse coding for denoising. In the Targeted Image Denoising (TID) method , targeted images are selected from a large dataset for each patch in the input noisy image for denoising, which is computationally expensive.

II-B Real-World Noisy Image Denoising

The methods emphasize much on the noise modeling, and they use Gaussian or MoG to model the noise in real-world noisy images. Nonetheless, the noise in real-world noisy images is very complex and hard to be modeled by explicit distributions . These works ignore the importance of learning image priors, which actually can be easier to model compared with modeling the complex realistic noise. In this paper, we propose a simple yet effective image prior learning method for real-world noisy image denoising. Due to its strong prior modeling ability, the proposed method simply models the noise as locally Gaussian, and it achieves highly competitive performance on real-world noisy image denoising.

III External Prior Guided Internal Prior Learning for Image Denoising

In this section, we first describe the learning of external prior, and then describe in detail the guided internal prior learning method, followed by the denoising algorithm.

Assume that a number of $L$ PGs are extracted from a set of external natural images, and the $l$-th PG is $\mathbf{\overline{X}}_{l}\triangleq\{\mathbf{\overline{x}}_{l,m}\}_{m=1}^{M},l=1,...,L$. A Gaussian Mixture Model (GMM) is learned to model the PG prior. The overall log-likelihood function is

The learning process is similar to the GMM learning in . Finally, a GMM model with $K$ Gaussian components is learned, and the learned parameters include mixture weights $\{\pi_{k}\}_{k=1}^{K}$, mean vectors $\{\bm{\mu}_{k}\}_{k=1}^{K}$, and covariance matrices $\{\bm{\Sigma}_{k}\}_{k=1}^{K}$. Note that the mean vector of each cluster is naturally zero, i.e., $\bm{\mu}_{k}=\bm{0}$.

To better describe the subspace of each Gaussian component, we perform singular value decomposition (SVD) on the covariance matrix:

The eigenvector matrices $\{\bm{U}_{k}\}_{k=1}^{K}$ will be employed as the external orthogonal dictionary to guide the internal sub-dictionary learning in next sub-section. The singular values in $\bm{S}_{k}$ reflect the significance of the singular vectors in $\bm{U}_{k}$. They will also be utilized as prior weights for weighted sparse coding in our denoising algorithm.

III-B Guided Internal Prior Learning

After the external PG prior model is learned from external natural clean images, we employ it to guide the internal PG prior learning for a given real-world noisy image. The guidance lies in two aspects. First, the external prior will guide the subspace clustering of internal noisy PGs. Second, the external prior will guide the orthogonal dictionary learning of internal noisy PGs.

Given a real-world noisy image $\mathbf{y}$, we extract $N$ (overlapped) local patches from it. Similar to the external prior learning stage, for the $n$-th ($n=1,...,N$) local patch we search its $M$ most similar (by Euclidean distance) patches around it to form a noisy PG, denoted by $\bm{Y}_{n}=\{\mathbf{y}_{n,1},...,\mathbf{y}_{n,M}\}$. Then the group mean of $\bm{Y}_{n}$, denoted by $\bm{\mu}_{n}$, is subtracted from each patch by $\bm{\overline{y}}_{n,m}\triangleq\mathbf{y}_{n,m}-\bm{\mu}_{n}$, leading to the mean subtracted noisy PG $\bm{\overline{Y}}_{n}\triangleq\{\bm{\overline{y}}_{n,m}\}_{m=1}^{M}$.

The external GMM prior models $\{\mathcal{N}(\bm{0},\bm{\Sigma}_{k})\}_{k=1}^{K}$ basically characterize the subspaces of natural high quality PGs. Therefore, we project each noisy PG $\bm{\overline{Y}}_{n}$ into the subspaces of $\{\mathcal{N}(\bm{0},\bm{\Sigma}_{k})\}_{k=1}^{K}$ and assign it to the most suitable subspace based on the posterior probability:

for $k=1,...,K$. Then $\bm{\overline{Y}}_{n}$ is assigned to the subspace with the maximum a-posteriori (MAP) probability $\max_{k}P(k|\bm{\overline{Y}}_{n})$.

III-B2 Guided Orthogonal Dictionary Learning

Assume that we have assigned all the internal noisy PGs $\{\bm{\overline{Y}}_{n}\}_{n=1}^{N}$ to their corresponding most suitable subspaces in $\{\mathcal{N}(\bm{0},\bm{\Sigma}_{k})\}_{k=1}^{K}$. For the $k$-th subspace, the noisy PGs assigned to it are $\{\bm{\overline{Y}}_{k_{n}}\}_{n=1}^{N_{k}}$, where $\bm{\overline{Y}}_{k_{n}}=[\bm{\overline{y}}_{k_{n},1},...,\bm{\overline{y}}_{k_{n},M}]$ and $\sum_{k=1}^{K}N_{k}=N$. We propose to learn an orthogonal dictionary $\bm{D}_{k}$ from each set of PGs $\bm{\overline{Y}}_{k_{n}}$ to characterize the internal PG prior with the guidance of the corresponding external orthogonal dictionary $\bm{U}_{k}$ (Eq. (2)). The reasons that we learn orthogonal dictionaries are two-fold. Firstly, the PGs $\{\bm{\overline{Y}}_{k_{n}}\}_{n=1}^{N_{k}}$ are in a subspace of the whole space of all PGs; therefore, there is no necessary to learn a redundant over-complete dictionary to characterize it, while an orthonormal dictionary has naturally zero mutual incoherence . Secondly, the orthogonality of dictionary can make the patch encoding in the testing stage very efficient, leading to an efficient denoising algorithm (please refer to sub-section III-C for more details).

We let the orthogonal dictionary $\bm{D}_{k}$ be

For notation simplicity, in the following development we ignore the subspace index $k$ for $\bm{\overline{Y}}_{k_{n}}$ and $\bm{D}_{k}$, etc. The learning of hybrid orthogonal dictionary $\bm{D}$ is performed under the following weighted sparse coding framework:

where $\bm{I}$ is the $3p^{2}$ dimensional identity matrix, $\bm{\alpha}_{n,m}$ is the sparse coding vector of the $m$-th patch $\bm{\overline{y}}_{n,m}$ in the $n$-th PG $\bm{\overline{Y}}_{n}$ and $\bm{\alpha}_{n,m,j}$ is the $j$-th element of $\bm{\alpha}_{n,m}$. $\lambda_{j}$ is the $j$-th regularization parameter defined as

where $\bm{S}_{k}(j)$ is the $j$-th singular value of diagonal singular value matrix $\bm{S}_{k}$ (please refer to Eq. (2)) and $\varepsilon$ is a small positive number to avoid zero denominator. Note that $\bm{D}_{\text{E}}=\bm{U}_{k}$ if $r=3p^{2}$ and $\bm{D}_{\text{E}}=\varnothing$ if $r=0$.

Updating Sparse Coding Coefficients: Given the orthogonal dictionary $\bm{D}^{(t)}$, we update each sparse coding vector $\bm{\alpha}_{n,m}$ by solving

Since dictionary $\bm{D}^{(t)}$ is orthogonal, the problems (7) has a closed-form solution

where $\bm{\lambda}=\frac{1}{2}[\lambda_{1},\lambda_{2},...,\lambda_{3p^{2}}]^{\top}$ is the vector of regularization parameter, $\text{sgn}(\bullet)$ is the sign function and $\odot$ means element-wise multiplication. The detailed derivation of Eq. (8) can be found in Appendix A.

Updating Internal Sub-dictionary: Given the sparse coding vectors $\{\bm{\alpha}_{n,m}^{(t+1)}\}$, we update the internal sub-dictionary by solving

The orthogonality of internal sub-dictionary $\bm{D}_{\text{I}}^{(t+1)}$ can be checked by $(\bm{D}_{\text{I}}^{(t+1)})^{\top}(\bm{D}_{\text{I}}^{(t+1)})=\bm{V}_{\text{I}}\bm{U}_{\text{I}}^{\top}\bm{U}_{\text{I}}\bm{V}_{\text{I}}^{\top}=\bm{I}_{(3p^{2}-r)}$. In fact, the Theorem 1 provides a sufficient and necessary condition to guarantee the existence of the closed-form solution for the internal sub-dictionary of the problem (9).

Besides, if $\text{rank}(\Sigma)=3p^{2}-r$, $\hat{\mathcal{D}}=\mathcal{U}\mathcal{V}^{\top}$ is also the sufficient condition of problem (11).

The proof of the Theorem 1 can be found in Appendix B. Though the problem (9) has a closed-form solution by SVD , the uniqueness of solution cannot be guaranteed since the matrices $(\bm{I}_{3p^{2}\times 3p^{2}}-\mathcal{E}\mathcal{E}^{\top})\mathcal{Y}\mathcal{A}^{\top}$ as well as $\mathcal{U}$ and $\mathcal{V}$ may be reduced to matrices of lower rank. Hence, we also analyze the uniqueness of the solution $\hat{\mathcal{D}}$ by the following Theorem 2, whose proof can be found in Appendix C.

The above alternative updating steps are repeated until the number of iterations exceeds a preset threshold. In each step, the energy value of the objective function (5) is decreased and we empirically found that the proposed model usually converges in 10 iterations. We summarize the procedures in Algorithm 1.

III-C The Denoising Algorithm

The denoising of the given noisy image $\mathbf{y}$ can be simultaneously done with the guided internal sub-dictionary learning process. Once we obtain the solutions of sparse coding vectors $\{\hat{\bm{\alpha}}_{n,m}^{(T)}\}$ in Eq. (8) and the orthogonal dictionary $\bm{D}^{(T)}=[\bm{D}_{\text{E}}\ \bm{D}_{\text{I}}^{(T)}]$ in Eq. (9), the latent clean patch $\hat{\mathbf{y}}_{n,m}$ of the $m$-th noisy patch in PG $\bm{Y}_{n}$ is reconstructed as

where $\bm{\mu}_{n}$ is the group mean of $\bm{Y}_{n}$. The latent clean image is then reconstructed by aggregating all the reconstructed patches in all PGs. We perform the above denoising procedures for several iterations for better denoising outputs. The proposed denoising algorithm is summarized in Algorithm 2.

IV Experiments

The noise in real-world images is very complex due to the many factors such as sensors, lighting conditions and camera settings. It is difficult to evaluate one algorithm by tuning its parameters for all these different settings. In this work, we fix the parameters of our algorithm and apply it to all the testing datasets, though they were captured by different types of sensors and under different camera settings. The parameters of our method include the patch size $p$, the number of similar patches $M$ in a patch group (PG), the window size $W$ for PG searching, the number of Gaussian components $K$ in GMM, the number of atoms $r$ in the external sub-dictionaries, the sparse regularization parameter $\lambda$, the iteration numbers $T$ for solving problem (5) and $IteNum$ for Alg. 2.

The performance of our proposed method varies little when we set patch size between $p=6$ and $p=9$, and we fix the patch size as $p=6$ to save computational cost. The search window is fixed to $W=31$ to balance computational cost and denoising accuracy of the proposed method. The number of patches in a patch group is set as $M=10$, while using more patches will not bring clear benefits. We learn the external GMM prior with 3.6 million PGs extracted from the Kodak PhotoCD Dataset (http://r0k.us/graphics/kodak/), which includes 24 high quality color images. The number of Gaussians in GMM is set as $K=32$, while using more Gaussians can only bring slightly better performance but cost more computational resources. The number of atoms in the external sub-dictionaries affects little the performance when it is set between $r=27$ and $r=81$, and we set it as $r=54$ to make the external and internal sub-dictionaries have the same number of atoms. We set the number of iterations as $T=2$ for solving the problem (5), while the number of iterations for Alg. 2 is set as $IteNum=4$.

One key parameter of our model is the regularization parameter $\lambda$. Fig. 3 plots the curves of PSNR/SSIM results w.r.t $\lambda$ on the 15 cropped image in dataset . One can see that our proposed method achieves good PSNR/SSIM performance within a certain range of $\lambda$. Similar observations can be made on other datasets. We fix $\lambda=0.001$ in the paper, and it works well across the three datasets used in our experiments.

All the parameters of our method are fixed in all experiments, which are run under the Matlab2014b environment on a machine with Intel(R) Core(TM) i7-5930K CPU of 3.5GHz and 32GB RAM. We will release the code with the publication of this work.

IV-B The Testing Datasets

We evaluate the proposed method on three real-world noisy image datasets, where the images were captured under indoor or outdoor lighting conditions by different types of cameras and camera settings.

Dataset 1. The first dataset is provided in , which includes noisy images of 11 static scenes. The noisy images were collected under controlled indoor environment. Each scene was shot 500 times under the same camera and camera setting. The mean image of the 500 shots is roughly taken as the “ground truth”, with which the PSNR and SSIM can be computed.

Since the image size is very large (about $7000\times 5000$) and the 11 scenes share repetitive contents, the authors of cropped 15 smaller images (of size $512\times 512$) to perform experiments. In order to evaluate the proposed methods more comprehensively, we cropped 60 images of size $500\times 500$ from the dataset for experiments. Some samples are shown in Fig. 4. Note that our cropped 60 images and the 15 cropped images by the authors of are from different shots.

Dataset 2 is called the Darmstadt Noise Dataset (DND) , which includes 50 different pairs of images of the same scenes captured by Sony A7R, Olympus E-M10, Sony RX100 IV, and Huawei Nexus 6P. The real-world noisy images are collected under higher ISO values with shorter exposure time, while the “ground truth” images are captured under lower ISO values with longer exposure times. Since the captured images are of megapixel-size, the authors cropped 20 bounding boxes of $512\times 512$ pixels from each image in the dataset, yielding $50\times 20=1000$ test crops in total. Some samples are shown in Fig. 5. Note that the “ground truth” images of this dataset have not been released yet, but one can submit the denoised images to the Project Website and get the average PSNR (dB) and SSIM results.

Dataset 3. On one hand, the scenes of Dataset 1 are mostly printed photos, and they cannot represent realistic objects and scenes with different reflectance properties. On the other hand, the Dataset 2 contains repetitive contents in the 20 cropped images for each of the 50 scenes. To remedy the limitations of Dataset 1 and Dataset 2, we construct another dataset which contains images of 10 different scenes captured by Canon 80D and Sony A7II cameras with more ISO settings and more comprehensive scenes. The ISO settings in our dataset are 800, 1600, 3200, 6400, 12800 while those of Dataset 1 are 1600, 3200, 6400. Compared to Dataset 2, our new dataset is more comprehensive on scene contents. Similar to Dataset 1, each scene was captured 500 shots, and the mean image of these 500 shots can be used a kind of ground-truth to evaluate the denoising algorithms. Fig. 6 shows some cropped images of the scenes in our dataset. One can see that the images contain a lot of different realistic objects with varying colors, shapes, materials, etc.

Our dataset provides real-world noisy images of realistic objects with different ISO settings. It can be used to more fairly evaluate the performance of different real-world noisy image denoising methods. Consider that the image resolution is very high (about $4000\times 4000$), for the convenience of experimental studies, we cropped 100 (10 for each scene) smaller images (of size $512\times 512$) from it to perform experiments. The whole dataset will be made publically available with the publication of this paper.

IV-C Comparison among external, internal and guided internal priors

To demonstrate the advantages of external prior guided internal prior learning, we perform real-world noisy image denoising by using external priors only (denoted by “External”), internal priors only (denoted by “Internal”), and the proposed guided internal priors (denoted by “Guided Internal”), respectively. For the “External” method, we utilize the full external dictionaries (i.e., $r=108$ in Eq. (5)) for denoising. For the “Internal” method, the overall framework is similar to the method of . A GMM model (with $K=32$ Gaussians) is directly learned from the PGs extracted from the given noisy image without using any external data, and then the internal orthogonal dictionaries are obtained via Eq. (2) to perform denoising. All parameters of the “External” and “Internal” methods are tuned to achieve their best performance.

We compare the three methods on the 60 cropped images from Dataset 1 . The average PSNR and run time are listed in Table I. The best results are highlighted in bold. It can be seen that “Guided Internal” method achieves better PSNR than both “External” and “Internal” methods. In addition, the “Internal” method is very slow because it involves online GMM learning, while the “Guided Internal” method is only a little slower than the “External” method. Figs. 7 and 8 show the denoised images of two noisy images by the three methods. One can see that the “External” method is good at recovering large-scale structures (see Fig. 7) while the “Internal” method is good at recovering fine-scale textures (see Fig. 8). By utilizing external priors to guide the internal prior learning, our proposed method can effectively recover both the large-scale structures and fine-scale textures.

IV-D Comparison with State-of-the-Art Denoising Methods

Comparison methods. We compare the proposed method with state-of-the-art image denoising methods, including GAT-BM3D , CBM3D , WNNM , TID , MLP , DnCNN , CSF , TNRD , Noise Clinic (NC) , Cross-Channel (CC) , and Neat Image (NI) . Among these methods, GAT-BM3D is a state-of-the-art Poisson noise reduction method. The method CBM3D is a state-of-the-art method for color image denoising and the noise on color images is assumed to be additive white Gaussian. The methods of WNNM, MLP, DnCNN, CSF, and TNRD are state-of-the-art Gaussian noise removal methods for grayscale images, and we apply them to each channel of color images for denoising. NC is a blind image denoising method, and NI is a set of commercial software for image denoising, which has been embedded into Photoshop and Corel PaintShop. The code of CC is not released but its results on the 15 cropped images are available at . Therefore, we only compare with it on the 15 cropped images in Dataset 1 .

Noise level of comparison methods. For the CBM3D method, the standard deviation of noise on color images should be given as a parameter. For methods of WNNM, MLP, CSF, and TNRD, the noise level in each color channel should be input. For the DnCNN method, it is trained to deal with noise in a range of levels $0\sim 55$. We retrain the models of discriminative denoising methods MLP, CSF, and TNRD (using the released codes by the authors) at different noise levels from $\sigma=5$ to $\sigma=50$ with a gap of $5$. The denoising is performed by processing each channel with the model trained at the same (or nearest) noise level. The noise levels ($\sigma_{r},\sigma_{g},\sigma_{b}$) in R, G ,B channels are assumed to be Gaussian and can be estimated via some noise estimation methods . In this paper, we employ the method to estimate the noise level for each color channel.

Results on Dataset 1. As described in section 4.2, there is a mean image for each of the 11 scenes used in Dataset 1 , and those mean images can be roughly taken as “ground truth” images for quantitative evaluation of denoising algorithms. We firstly perform quantitative comparison on the 15 cropped images used in . The results on PSNR (dB) and speed (second) of GAT-BM3D, CBM3D, WNNM, TID, MLP, CSF, TNRD, DnCNN, NC, NI and CC are listed in Table II (The results of CC are copied from the original paper ). The best PSNR results of each image are highlighted in bold. One can see that on 8 out of the 15 images, our method achieves the best PSNR values. CC achieves the best PSNR on 3 of the 15 images. It should be noted that in the CC method, a specific model is trained for each camera and camera setting, while our method uses the same model for all images. On average, our proposed method has 0.27dB PSNR improvements over the second best method CC and much higher PSNR gains over other competing methods. The method GAT-BM3D does not work well on most images. This is because real world noise is much more complex than Poisson.

Figs. 9 and 10 show the denoised images of one scene captured by Canon 5D Mark 3 at ISO = 3200 and Nikon D800 at ISO = 6400, respectively. We can see that GAT-BM3D, CBM3D, TID, DnCNN, NC, NI and CC would either remain noise or generate artifacts, while TNRD over-smooths much the image. By using the external prior guided internal priors, our proposed method preserves edges and textures better than other methods while removing the noise, leading to visually more pleasant outputs. Specifically, Fig. 10 is used to illustrate the denoising performance of our method on fine-scale textures such as hair, which is a very challenging task. Even the “ground truth” mean image cannot show very clear details of the hair. Though our method cannot reproduce clearly the details (e.g., the local direction of hair in some regions), it demonstrates the best visual results among the competing methods. More comparisons on visual quality and SSIM index can be found in the supplementary file.

We then perform denoising experiments on the 60 images we cropped from . The average PSNR results are listed in Table III (CC is not compared since the code is not available). Again, our proposed method achieves much better PSNR results than the other methods. The improvements of our method over the second best method (TNRD) are 0.43dB on PSNR. Fig. 11 shows the denoised images of one scene captured by Nikon D800 at ISO = 3200. We can see again that the proposed method obtain better visual quality than other competing methods. More comparisons on visual quality and SSIM can be found in the supplementary file.

Results on Dataset 2. In Table IV, we list the average PSNR (dB) results of the competing methods on the 1000 cropped images in the DND dataset . We can see again that the proposed method achieves better performance than the other competing methods. Note that the “ground truth” images of this dataset have not been released yet, so we are not able to calculate the PSNR and SSIM results for each noisy image in this dataset, nor compare with the “ground truth” mean image. However, one can submit the denoised images to the project website and get the average PSNR and SSIM results on the whole 1000 images. Fig. 12 shows the denoised images of a scene “0001_2” captured by a Nexus 6P phone . The noise level in this image is relatively high. Hence, this image can be used to justify the performance of the proposed method on real-world noisy images with lower PSNR (around 20dB). One can see that the proposed method achieves visually more pleasing results than the other denoising methods. More comparisons on visual quality and SSIM can be found in the supplementary file.

Results on Dataset 3. Similar to Dataset 1 , there is a “ground truth” image for each of the 10 scenes used in our constructed Dataset 3. We perform quantitative comparison on the 100 cropped images. The average PSNR results of competing methods are listed in Table IV. We can see that our proposed method achieves much better PSNR results than the other methods. The improvements of our method over the second best method (TNRD) is 0.16dB on PSNR. Fig. 13 shows the denoised images of one scene captured by Canon 80D at ISO = 12800. We can see again that the proposed method removes the noise while maintains better details (such as the vertical black shadow area) than other competing methods. More comparisons on visual quality and SSIM can be found in the supplementary file.

Comparison on speed. Efficiency is an important aspect to evaluate the efficiency of algorithms. We compare the speed of all competing methods except for CC. All experiments are run under the Matlab2014b environment on a machine with Intel(R) Core(TM) i7-5930K CPU of 3.5GHz and 32GB RAM. The average running time (second) of the compared methods on the 100 real-world noisy images is shown in Table V. The least average running time are highlighted in bold. One can easily see that the commercial software Neat Image (NI) is the fastest method with highly optimized code. For a $512\times 512$ image, NI costs about 0.6 second. The other methods cost from 5.2 (TNRD) to 152.2 (WNNM) seconds, while the proposed method costs about 24.1 seconds. It should be noted that GAT-BM3D, CBM3D, TNRD, and NC are implemented with compiled C++ mex-function and with parallelization, while WNNM, TID, MLP, CSF, DnCNN, and the proposed method are implemented purely in Matlab.

V Conclusion

We proposed a new prior learning method for the real-world noisy image denoising problem by exploiting the useful information in both external and internal data. We first learned Gaussian Mixture Models (GMMs) from a set of clean external images as general image prior, and then employed the learned GMM model to guide the learning of adaptive internal prior from the given noisy image. Finally, a set of orthogonal dictionaries were output as the external-internal hybrid prior models for image denoising. Extensive experiments on three real-world noisy image datasets, including a new dataset constructed by us by different types of cameras and camera settings, demonstrated that our proposed method achieves much better performance than state-of-the-art image denoising methods in terms of both quantitative measure and visual perceptual quality.

Appendix A Closed-Form Solution of the Weighted Sparse Coding Problem (7)

For notation simplicity, we ignore the indices $n,m,t$ in problem (7). It turns into the following weighted sparse coding problem:

Since $\bm{D}$ is an orthogonal matrix, problem (14) is equivalent to:

For simplicity, we denote $\mathbf{z}=\mathbf{D^{T}y}$. Here we have $\lambda_{j}>0$, $j=1,...,3p^{2}$, then problem (15) can be written as:

The problem (16) is separable w.r.t. each $\bm{\alpha}_{j}$ and hence can be simplified to $3p^{2}$ independent scalar minimization problems:

where $j=1,...,3p^{2}$. Taking derivative of $\bm{\alpha}_{j}$ in problem (17) and setting the derivative to be zero. There are two cases for the solution.

(a) If $\bm{\alpha}_{j}\geq 0$, we have $2(\bm{\alpha}_{j}-\mathbf{z}_{j})+\lambda_{j}=0,$ and the solution is $\hat{\bm{\alpha}}_{j}=\mathbf{z}_{j}-\frac{\lambda_{j}}{2}\geq 0.$ So $\mathbf{z}_{j}\geq\frac{\lambda_{j}}{2}>0$, and the solution $\hat{\bm{\alpha}}_{j}$ can be written as $\hat{\bm{\alpha}}_{j}=\text{sgn}(\mathbf{z}_{j})*(|\mathbf{z}_{j}|-\frac{\lambda_{j}}{2}),$ where $\text{sgn}(\bullet)$ is the sign function.

(b) If $\bm{\alpha}_{j}<0$, we have $2(\bm{\alpha}_{j}-\mathbf{z}_{j})-\lambda_{j}=0$ and the solution is $\hat{\bm{\alpha}}_{j}=\mathbf{z}_{j}+\frac{\lambda_{j}}{2}<0.$ So $\mathbf{z}_{j}<-\frac{\lambda_{j}}{2}<0$, and the solution $\hat{\bm{\alpha}}_{j}$ can be written as $\hat{\bm{\alpha}}_{j}=\text{sgn}(\mathbf{z}_{j})*(-\mathbf{z}_{j}-\frac{\lambda_{j}}{2})=\text{sgn}(\mathbf{z}_{j})*(|\mathbf{z}_{j}|-\frac{\lambda_{j}}{2}).$

In summary, we have the final solution of the weighted sparse coding problem (14) as:

where $\bm{\lambda}=\frac{1}{2}[\lambda_{1},\lambda_{2},...,\lambda_{3p^{2}}]^{\top}$ is the vector of regularization parameter and $\odot$ means element-wise multiplication.

Appendix B Proof of the Theorem 1

We firstly prove the necessary condition. Since $\mathcal{D}^{\top}\mathcal{D}=\bm{I}_{(3p^{2}-r)\times(3p^{2}-r)}$, we have

The Lagrange function is $\mathcal{L}=\text{Tr}(\mathcal{A}\mathcal{Y}^{\top}\mathcal{D})-\text{Tr}(\Gamma_{1}(\mathcal{D}^{\top}\mathcal{D}-\bm{I}_{(3p^{2}-r)\times(3p^{2}-r)}))-\text{Tr}(\Gamma_{2}(\mathcal{D}^{\top}\mathcal{E}))$, where $\Gamma_{1}$ and $\Gamma_{2}$ are the Lagrange multipliers. Take the derivative of $\mathcal{L}$ w.r.t. $\mathcal{D}$ and set it to be matrix $\bm{0}$ of conformal dimensions, we can get

Since $\mathcal{D}^{\top}\mathcal{D}=\bm{I}_{(3p^{2}-r)\times(3p^{2}-r)}$ and $\mathcal{E}^{\top}\mathcal{D}=\bm{0}_{3p^{2}\times(3p^{2}-r)}$, by left multiplying both sides of the Eq. (22) by $\mathcal{E}^{\top}$, we have

Put the Eq. (22) back into Eq. (21), we have

Right multiplying both sides of Eq. (23) by $\mathcal{D}^{\top}$, we have

This shows that $(\bm{I}_{3p^{2}\times 3p^{2}}-\mathcal{E}\mathcal{E}^{\top})\mathcal{Y}\mathcal{A}^{\top}\mathcal{D}^{\top}$ is a symmetric matrix of order $3p^{2}\times 3p^{2}$. Then we perform economy (or reduced) singular value decomposition (SVD) on $(\bm{I}_{3p^{2}\times 3p^{2}}-\mathcal{E}\mathcal{E}^{\top})\mathcal{Y}\mathcal{A}^{\top}=\mathcal{U}\Sigma\mathcal{V}^{\top}$, there is

Hence, we have $\mathcal{U}=\mathcal{D}\mathcal{V}$, or equivalently $\hat{\mathcal{D}}=\mathcal{U}\mathcal{V}^{\top}$. The necessary condition is proved.

Now we prove the sufficient condition. If $\hat{\mathcal{D}}=\mathcal{U}\mathcal{V}^{\top}$, then $\hat{\mathcal{D}}^{\top}\hat{\mathcal{D}}=\bm{I}_{(3p^{2}-r)\times(3p^{2}-r)}$. To prove $\mathcal{E}^{\top}\hat{\mathcal{D}}=\bm{0}_{3p^{2}\times(3p^{2}-r)}$, we left multiply both sides of Eq. (25) by $\mathcal{E}^{\top}$ and have $\bm{0}_{3p^{2}\times(3p^{2}-r)}=\mathcal{E}^{\top}(\bm{I}_{3p^{2}\times 3p^{2}}-\mathcal{E}\mathcal{E}^{\top})\mathcal{Y}\mathcal{A}^{\top}\hat{\mathcal{D}}^{\top}=\mathcal{E}^{\top}\mathcal{U}\Sigma\mathcal{V}^{\top}\hat{\mathcal{D}}^{\top}=\mathcal{E}^{\top}\mathcal{U}\Sigma\mathcal{U}^{\top}$ . It means that $\mathcal{E}^{\top}\mathcal{U}\Sigma\mathcal{U}^{\top}=\bm{0}_{3p^{2}\times 3p^{2}}$. This only happens when $\mathcal{E}^{\top}\mathcal{U}=\bm{0}_{3p^{2}\times(3p^{2}-r)}$ since $\text{rank}(\Sigma)=3p^{2}-r$ and $\mathcal{U}\Sigma\mathcal{U}^{\top}$ is positive definite. Then $\mathcal{E}^{\top}\hat{\mathcal{D}}=\mathcal{E}^{\top}\mathcal{U}\mathcal{V}^{\top}=\bm{0}_{3p^{2}\times(3p^{2}-r)}$.

Finally we prove that $\hat{\mathcal{D}}=\mathcal{U}\mathcal{V}^{\top}$ is the solution of

Note that by cyclic perturbation which retains the trace unchanged and due to $\mathcal{E}^{\top}\hat{\mathcal{D}}=\bm{0}_{3p^{2}\times(3p^{2}-r)}$, we have $\text{Tr}(\mathcal{Y}^{\top}\hat{\mathcal{D}}\mathcal{A})=\text{Tr}(\mathcal{Y}\mathcal{A}^{\top}\hat{\mathcal{D}}^{\top})=\text{Tr}((\bm{I}_{3p^{2}\times 3p^{2}}-\mathcal{E}\mathcal{E}^{\top})\mathcal{Y}\mathcal{A}^{\top}\hat{\mathcal{D}}^{\top})=\text{Tr}(\mathcal{U}\Sigma\mathcal{V}^{\top}\mathcal{V}\mathcal{U}^{\top})=\text{Tr}(\Sigma).$ For every $\mathcal{D}$ satisfying that $\mathcal{D}^{\top}\mathcal{D}=\bm{I}_{(3p^{2}-r)\times(3p^{2}-r)}$, $\mathcal{E}^{\top}\mathcal{D}=\bm{0}_{3p^{2}\times(3p^{2}-r)}$, we have $\text{Tr}(\mathcal{Y}^{\top}\mathcal{D}\mathcal{A})=\text{Tr}((\bm{I}_{3p^{2}\times 3p^{2}}-\mathcal{E}\mathcal{E}^{\top})\mathcal{Y}\mathcal{A}^{\top}\mathcal{D}^{\top})=\text{Tr}(\mathcal{U}\Sigma\mathcal{V}^{\top}\mathcal{D}^{\top})=\text{Tr}(\Sigma\mathcal{V}^{\top}\mathcal{D}^{\top}\mathcal{U})$. By using a generalization version of the Kristof’s Theorem , we have $\text{Tr}(\mathcal{Y}^{\top}\mathcal{D}\mathcal{A})=\text{Tr}(\Sigma\mathcal{V}^{\top}\mathcal{D}^{\top}\mathcal{U})\leq\text{Tr}(\Sigma).$ The equality is obtained at $\mathcal{V}^{\top}\mathcal{D}^{\top}\mathcal{U}=\bm{I}_{(3p^{2}-r)\times(3p^{2}-r)}$, i.e., $\mathcal{D}=\mathcal{U}\mathcal{V}^{\top}=\hat{\mathcal{D}}$. This completes the proof. ∎

Appendix C Proof of the Theorem 2

Before we prove the Theorem 2, we need firstly prove the following Lemma 1.

Since $\text{rank}(\mathcal{E}\mathcal{E}^{\top})\leq\min\{\text{rank}(\mathcal{E}),\text{rank}(\mathcal{E}^{\top})\}=r$ and $\text{rank}(\mathcal{E}\mathcal{E}^{\top})\geq\text{rank}(\mathcal{E})+\text{rank}(\mathcal{E}^{\top})-r=2r-r=r$ by Sylvester’s inequality, we have $\text{rank}(\mathcal{E}\mathcal{E}^{\top})=r$. Then, $\text{rank}(\bm{I}_{3p^{2}\times 3p^{2}}-\mathcal{E}\mathcal{E}^{\top})\geq\text{rank}(\bm{I}_{3p^{2}\times 3p^{2}})-\text{rank}(\mathcal{E}\mathcal{E}^{\top})\geq 3p^{2}-r$. ∎

The $\text{rank}(\Sigma)$ ($\Sigma$ is defined in Theorem 1) depends on $\text{rank}(\bm{I}_{3p^{2}\times 3p^{2}}-\mathcal{E}\mathcal{E}^{\top})$, $\text{rank}(\mathcal{Y})$ and $\text{rank}(\mathcal{A})$. Note that $\text{rank}(\mathcal{Y})\geq M$ and $\text{rank}(\mathcal{A})\geq\min\{3p^{2},M\}$ and $\text{rank}(\bm{I}_{3p^{2}\times 3p^{2}}-\mathcal{E}\mathcal{E}^{\top})\geq 3p^{2}-r$. Hence, $\text{rank}(\Sigma)\leq\min\{3p^{2}-r,M\}$.

Right multiplying both sides of Eq. (28) by $\mathcal{D}\mathcal{V}$ and left multiplying each side by $\mathcal{U}^{\top}$, we have

Thus, we have $\mathcal{D}=\mathcal{U}\Delta\mathcal{V}^{\top}$. That is, if $\text{rank}(\Sigma)<3p^{2}-r$, once we get the solution of $\hat{\mathcal{D}}=\mathcal{U}\mathcal{V}^{\top}$ in problem (19), $\mathcal{D}=\mathcal{U}\Delta\mathcal{V}^{\top}$ with suitable $\Delta$ is also the solution of problem (19). In fact, the number of solutions $\hat{\mathcal{D}}$ for problem (19) is $2^{3p^{2}-r-\text{rank}(\Sigma)}$ given fixed $\mathcal{U}$ and $\mathcal{V}$. ∎

References