Noise2Score: Tweedie’s Approach to Self-Supervised Image Denoising without Clean Images
Year: 2021

Abstract:
Recently, there has been extensive research interest in training deep networks to denoise images without clean reference. However, the representative approaches such as Noise2Noise, Noise2Void, Stein’s unbiased risk estimator (SURE), etc. seem to differ from one another and it is difficult to find the coherent mathematical structure. To address this, here we present a novel approach, called Noise2Score, which reveals a missing link in order to unite these seemingly different approaches. Specifically, we show that image denoising problems without clean images can be addressed by finding the mode of the posterior distribution and that the Tweedie’s formula offers an explicit solution through the score function (i.e. the gradient of loglikelihood). Our method then uses the recent finding that the score function can be stably estimated from the noisy images using the amortized residual denoising autoencoder, the method of which is closely related to Noise2Noise or Nose2Void. Our Noise2Score approach is so universal that the same network training can be used to remove noises from images that are corrupted by any exponential family distributions and noise parameters. Using extensive experiments with Gaussian, Poisson, and Gamma noises, we show that Noise2Score significantly outperforms the state-of-the-art self-supervised denoising methods in the benchmark data set such as (C)BSD68, Set12, and Kodak, etc.

heading: 1 Introduction
text: Bayesian inference, which derives the posterior probability using a prior probability and a likelihood function for the observed data, has been an important tool in statistics. This approach has been used extensively for image denoising from early ages to the modern era of deep learning. For example, in the recent unsupervised deep learning approach for image denoising using Stein’s risk estimate (SURE) [1], the unknown Bayesian risk is replaced by the SURE that can be calculated from the noisy measurement so that deep neural network training is performed by minimizing it. Unfortunately, this method is sensitive to hyper-parameters, and the neural network must be retrained if the underlying noise model varies [2]. On the other hand, there has been increased research interest in image denoisers that can be trained by minimizing variants of empirical risks that are not associated with clean data. Noise2Noise [3] was the first representative approach that does not require clean data. Unfortunately, multiple noisy versions of the same images are necessary for training. To address this, self-supervised learning approaches such as Noise2Void [4], Noise2Self [5], etc. have been developed in order to use only a single noisy image. This class of approaches, which we will call Noise2X throughout the paper, are especially important for practical applications, where noiseless clean images or multiple noisy realization of the same image are difficult or impossible to collect.
35th Conference on Neural Information Processing Systems (NeurIPS 2021).
Though SURE and Noise2X hold promise for practical applications, one of the fundamental questions is how these seemingly different approaches are related and whether there is a coherent mathematical theory that can be leveraged to further improve the performance. Although we are aware of a preliminary prior work [6] attempting to unite a subset of them, there is largely a lack of a principled method of designing image denoiser without clean data.
In this article, we try to approach this open problem from Bayesian statistics in a completely different route. Instead of minimizing different forms of empirical risks or SURE, one of the most important contributions of this paper is the discovery of the importance of classical result from Bayesian statistics - the Tweedie’s formula [7] which provide an explicit way of computing the posterior mean of canonical parameters from the noisy measurements corrupted with exponential family noises. Specifically, we show that the Tweedie’s formula provides a unified approach for image denoising from any exponential family noises through a score function (i.e. the gradient of the loglikelihood). Therefore, the self-supervised image denoising problem without clean image can be reduced to the problem of estimating the score function.
This change of perspective has many important theoretical implications as well as the flexibility with regard to the implementation of the algorithm. In fact, the score function estimation problem has been an important research topic in Bayesian statistics and machine learning [8–10]. In particular, Alain and Bengio [10] showed that the minimization of the denoising autoencoder (DAE) objective function provides an explicit way of approximating the score function. This result was further extended using the amortized residual denoising autoencoder (AR-DAE) for numerical stability and accuracy [11]. Interestingly, the training method for DAE or its variants is strikingly similar to that of Noise2X. Therefore, by combining this with the Tweedie’s formula, we can obtain a novel unified framework of Bayesian approach for self-supervised image denoising. The conceptual difference and similarity of our method to the existing approaches are illustrated in Fig. 1(a)(b). As the first step of our algorithm is similar to Noise2X, we call our method Noise2Score.
Noise2Score is so powerful that it can be used to deal with any exponential family noises. Moreover, in contrast to the SURE approach, in which the network must be retrained using a different loss function depending on the noise model and parameters [1, 12], Noise2Score is universal in the sense that the identical neural network training step is used regardless of noise models. This property can be exploited to extend Noise2Score to blind setups in which noise parameters are unknown and should be estimated with minimal complexity. In addition to the novel theoretical findings mentioned above, our empirical results using additive Gaussian, Poisson and Gamma noise models have shown that Noise2Score significantly outperforms SURE and Noise2X under similar experimental conditions.

heading: 2 Related Works
text: Here we give a brief overview of the existing deep learning approaches to denoise images. As shown below, their main differences arise from the different choices of the distance metric d(·, ·) in Fig. 1.

heading: 2.1 Image denoising using supervised learning
text: Let X and Y denotes the spaces of clean and noisy images with the probability measure PX and PY , respectively. When the paired ground-truth data x is available for a noisy measurement y, i.e.
x := x(y), then the supervised training can be carried out by minimizing the following loss function:
`sup(Θ) = Ey∼PY d(x, FΘ(y)), where d(x, FΘ(y)) := ‖x− FΘ(y)‖2 (1) where FΘ(y) is a deep neural network parameterized by Θ with input y and here we use l2 distance metric to simplify the explanation. Although this approach has significantly improved performance compared to the classical denoising algorithms [13], the supervised learning approach has a drawback as the large number of paired clean images are necessary. In real world applications, the acquisition for paired clean images are difficult or even impossible.

heading: 2.2 Noise2X for image denoising without clean images
text: To address this issue, various forms of the denoisers that can be trained without clean images have been proposed, which includes Noise2Noise [3], Noise2Void [4], Noise2Self [14], Noise2Same [5], etc. These approaches are usually based on the variants of loss functions that do not depend on clean images:
`Noise2X(Θ) = Ey∼PY d(y′, FΘ(y)), where d(y′, FΘ(y)) := ‖y′ − FΘ(y)‖2 (2) where the target image y′ is related to y in unique ways depending on algorithms. For example, in Noise2Noise [3], y′ is another noisy realization of the same underlying clean images. However, Noise2Noise is not practical since multiple noisy realizations of the same image are not usually available in practice. Noise2Void [14] generates y′ by adding blind spots on y. In Noise2Self [14], so-called J-invariant mask J is added on the input y, and y′ is defined as the image on the complementary mask. All of these changes of the target images are intended to prevent the network from converging to a trivial identity mapping.

heading: 2.3 Image denoising using Stein’s Unbiased Risk Estimate (SURE)
text: Rather than using the variants of a loss function as in Noise2X, Soltanayev et al [1] proposed a training scheme using the loss function from Stein’s unbiased risk estimate (SURE) [2] which has additional regularization term. Specifically, the loss function for the SURE denoiser is defined by
`SURE(Θ) = Ey∼PY d(y, FΘ(y)), d(y, FΘ(y)) := ‖y − FΘ(y)‖2 + 2σ2divyFΘ(y) (3) Here, the second term in d(y, FΘ(y)) of (3) is the regularization term to prevent the network from converging to a trivial identity mapping, where divy denotes the divergence with respect to y. Due to the difficulty of calculating the divergence term, the authors in [1] relied on MonteCarlo SURE [15] which calculates the divergence term using MonteCarlo simulation. This introduces additional hyperparameters, on which the final results critically depend. Yet another limitation of SURE image denoiser is that the loss function is determined by the noise model. For example, in [12], the authors derived Poisson Unbiased risk estimator (PURE) to train deep CNN for dealing with Poisson noise, whose loss function is completely different from (3). Even under the same Gaussian noise models, if the noise variance σ2 varies, the neural network FΘ must be trained again, which causes additional computational burden. This differs from Noise2X, which uses the same loss function regardless of the noise models.

heading: 3 Main Contribution: Noise2Score
text: In contrast to the existing approaches shown in Fig. 1(a), our denoising approach is composed of two steps as illustrated in Fig. 1(b), which is inspired by the results from Bayesian statistics. In the following, we explain the details of each step.

heading: 3.1 Tweedie’s formula of the posterior mean for image denoising
text: Suppose that the noisy measurement y is given by
y = x+ w, where w ∼ N (0, σ2I) (4) where x is the clean image, and N (0, σ2I) denotes the independent and identically distributed (i.i.d) Gaussian distribution with zero mean and variance σ2. Then, Tweedie’s formula calculates the posterior expectation of x given y as [16]
E[x|y] = y + σ2l′(y), where l′(y) = ∇y log p(y) (5)
where p(y) denotes the marginal distribution of y. For the image corrupted by Gaussian noises, the posterior mean is the minimum mean square error (MMSE) estimator, so that the denoised image can be obtained using (5), as long as we know l′(y) - the score function [17].
Tweedie’s formula was extended for general class of exponential family distribution [7]. Specifically, exponential family of probability distributions are defined as those distributions whose density have the following general form:
p(y|η) = p0(y) exp ( η>T (y)− ϕ(η) ) (6)
where the superscript > denotes the transpose operation. Here η is a canonical (vector) parameter of the family, T (y) is a (vector) function of y, ϕ(η) the cumulant generating function which makes p(y|η) integrate to 1, and p0(y) the density up to a scale factor when η = 0. Exponential family distribution covers a large class of important distributions such as the Gaussian, binomial, multinomial, Poisson, gamma, and beta distributions, as well as many others. Table 1 summarizes the corresponding parameters for each exponential family distribution.
Then, using the Bayes’ rule, the posterior density of η given y is given by [7]:
p(η|y) = exp ( η>T (y)− λ(y) ) [ p(η)e−ϕ(η) ] , where λ(y) = log ( p(y)
p0(y)
) (7)
where p(y) and p(η) denote the marginal distribution of y and η, respectively. This implies that the posterior density is again an exponential family distribution, and the mode of the posterior distribution can be obtained by finding the maximum of p(η|y). Specifically, by computing the gradient of log p(η|y) with respect to y and setting it to zero, the posterior estimate of the canonical parameter η̂ should satisfy the following equality:
η̂>T ′(y) = λ′(y) = −∇y log p0(y) +∇y log p(y) = −l′0(y) + l′(y) (8)
where l′(y) := ∇y log p(y) and l′0(y) := ∇y log p0(y) are score functions, and T ′(y) = ∇yT (y). In particular, if T (y) = y, we can obtain the following closed form solution for the posterior mean:
η̂ := E[η|y] = λ′(y) = −l′0(y) + l′(y) (9)
In the following, we derive denoising algorithms using Tweedie’s formulation for several representative distributions.
Additive Gaussian noise For Gaussian distribution, in the first row of Table 1, we have
p0(y) = 1√ 2πσ e− y2 2σ2 , η = x σ2 (10)
Accordingly,
E[η|y] = x̂ σ2 = y σ2 + l′(y), where x̂ := E[x|y] (11)
since l′0(y) = y/σ 2. Therefore, the posterior mean is given by
x̂ = y + σ2l′(y) (12)
which is equal to (5). Even with the different parameterization of Gaussian as shown in the second row of Table 1, we can arrive at the same result using (8).
Poisson noise For the case of Poisson noises, we have
p0(y) = 1
y! , l′0(y) = − ∇yΓ(y + 1) Γ(y + 1)
' − log ( y + 1
2
) (13)
where Γ(y) denotes the gamma function and the last approximation comes from [7]. Accordingly, we have
η̂ = log(x̂) = log ( y + 1
2
) + l′(y) (14)
This leads to the following posterior estimate of the image
x̂ = ( y + 1
2
) exp(l′(y)) (15)
This corresponds to Tweedie’s formula for Poisson case [16].
This result can be generalized to the low photon count sensing scenario, where the sensor measurement can be described by [18]: y = ζz, z ∼ Poisson(x/ζ) (16) where z is a random variable that follows the Poisson distribution in Table 1 and ζ > 0 is the gain of the acquisition process that is related to the noise level. Then, we have
p0(z) = 1
z! , η = log
( x
ζ
) (17)
Therefore, we have
x̂ ζ =
( z + 1
2
) exp(l′(z)) −→ x̂ = ( y + ζ
2
) el ′( yζ ) (18)
Gamma noise Gamma noise distribution can be used to model the speckle noises in various imaging application. Specifically, the image corrupted with the speckle noise is represented by
y = xn, n ∼ p(n;α, β) = β α
Γ(α) nα−1 exp(−βn) (19)
where p(n;α, β) denotes the Gamma distribution with (α, β) parameters. Here, (α, β) are the parameters which determine the noise level of Gamma noise distribution. For example, for the k-look measurement case, α = β = k [19]. Then, the probability p(y|x) can be obtained by replacing n in p(n;α, β) with n = y/x, and the resulting density function can be bound in Table 1. Then, we can easily show that
po(y) = 1, η = [ α− 1 −β/x ] , T (y) = [ log y y ] (20)
Using (8), we have
α− 1 y − β x = l′(y) =⇒ x̂ = βy (α− 1)− yl′(y)
(21)
Similar derivations for other exponential family distribution can be found in Table 1.

heading: 3.2 Score function estimation
text: So far, our derivation assumes the prior knowledge of the score function. In practice, this should be estimated, so here we describe how the score function can be estimated using a neural network.
Historically, Hyvärinen et al [8] was the first to derive a remarkable implicit score matching objective that no longer requires having an explicit score target but is nevertheless equivalent to the original problem that minimizes the expected quadratic distance between the model function and the score function of data. With the advance of the Denoising Auto Encoder (DAE) [20], the author of [9] observed that the minimization of the DAE objective function with the residual form is related to the score matching between the model and perturbed data. This observation was rigorously analyzed
by Alain and Bengio [10], who showed that as the perturbed noise becomes sufficiently small, DAE leads to the score function.
More specifically, in the denoising autoencoder (DAE), the following loss is minimized:
`DAE(Θ) = E y∼PY
u∼N (0,I) σa∼N (0,δ2)
‖y − FΘ(y + σau)‖2 (22)
According to [10], the optimal DAE FΘ∗(x) can be represented by
FΘ∗(y) = y + σ 2 al ′(y) + o(σ2a), (23)
where l′(y) is the score function defined in (5) and o(·) denotes the small “o” notation. In other words, for sufficiently small σa, we can approximate the score function as:
l̂′(y) = FΘ∗(y)− y
σ2a (24)
One of the downsides of the score function estimation using DAE is that the score function estimate in (24) is numerically unstable as σa → 0. To address the numerical instability and reduce the approximation error, Lim et al [11] recently proposed so-called the amortized residual DAE (ARDAE). Specifically, AR-DAE is trained by minimizing the following the objective function:
`AR−DAE(Θ) = E y∼PY
u∼N (0,I) σa∼N (0,δ2)
‖u+ σaRΘ(y + σau)‖2 (25)
where RΘ is from the residual form of the DAE:
FΘ(y) = σ 2 aRΘ(y) + y (26)
By plugging (26) in (25), we can obtain the original DAE cost in (22) up to a scale factor. Furthermore, using (24), we have
l̂′(y) = FΘ∗(y)− y
σ2a = RΘ(y) (27)
Therefore, the neural network trained with the AR-DAE (25) is a direct and stable way of estimating the score function. Therefore, we employ the AR-DAE to estimate the score function from the noisy measurement as the first step of our method (see Fig. 1(b)).

heading: 3.3 Relation to Noise2X and SURE
text: Although our Noise2Score was derived from a completely different perspective, it turns out that it has very important connections to Noise2X and SURE.
Specifically, if the noisy image y is corrupted by additive Gaussian, using Tweedie’s fomula for Gaussian noises in Table 1, the denoised image can be approximated by
x̂ = y + σ2 l̂′(y) = FΘ∗(y) (28)
where we use (27). Recall that FΘ∗(y) is the neural nework trained with DAE loss function in (22). The DAE training with (22) is basically inserting noises into the images and then find the mapping that removes the noises. In terms of adding extra noise to a noisy image, our method is closely related to Noisier2Noise [21]. However, in the training phase of Noisier2Noise, the main assumption is that noise from the same noise statistic as the original noisy image should be injected into the noisy images so that the neural networks learn to reduce the noise from the noisy input image. In our Noise2Score training method, however, Gaussian noise is added to noisy images at different noise levels regardless of the noise statistics to estimate the score function rather than noise. In addition, in contrast to Noisier2Noise, our method has a post-processing step with the Tweedie’s formula, which is determined by the noise model such as Poisson, Gamma, etc.
Yet another important connection is its relation to SURE. Using the residual representation in (26), the SURE cost function in (3) can be equivalently represented by
`SURE(Θ) = Ey∼PY { ‖y − FΘ(y)‖2 + 2σ2divyFΘ(y) } = Ey∼PY { σ4‖RΘ(y)‖2 + 2σ4divyRΘ(y) } + 2σ2dim(y) (29)
where dim(y) is the dimension of the vector y, which is constant. Up to the scaling factor, this cost function is identical to the implicit score matching objective by Hyvärinen et al [8], which is given by
`ISM (Θ) = Ey∼PY { 1
2 ‖ΨΘ(y)‖2 + divyΨΘ(y)
} (30)
where ΨΘ(y) is the score function estimate parameterized by Θ. Therefore, using the residual transform in (26), the SURE objective is nothing but to find the score function, which corresponds to the first step of our Noise2Score. Furthermore, due to the use of (26), SURE is not optimal for noise models other than Gaussian. Although there are variants of SURE for other noise type [22], the cost function of this method varies depending on the type of noise. However, as discussed in detail later, our method uses the same cost function for neural network training, but only the post-processing step differs depending on the type of noise. Moreover, considering that the modern research trend of the score function estimation has quickly evolved from the implicit score matching to DAE or AR-DAE due to the numerical instability and difficulty of computing the divergence term in (30), we can easily expect that the SURE denoiser also suffers from similar limitations.

heading: 3.4 Universal neural estimation for blind denosing problems
text: Recall that Noise2Score is composed of two steps: the score function estimation by minimizing the AR-DAE loss function, after which the final denoised image is obtained by applying Tweedie’s formula shown in Table 1. This decoupling allows Noise2Score to have unique advantages compared to the conventional approaches.
Specifically, one of the most important advantages of Noise2Score is that AR-DAE training is universal in the sense that the same loss function is used regardless of noise models and parameters. More specifically, in (25), the model and noise parameter dependency is only through the sampling y from PY . This sampling step can be performed by random patch cropping from the input images, which procedure is identical regardless of PY being Gaussian, Poisson, Gamma, etc.
Yet another advantage is that this property can be utilized to estimate the unknown noise parameters without retraining neural networks. For example, in the case of Gaussian denoising in Table 1, we have the following form of the denoised images:
x̂(σ) = y + σ2 l̂′(y) (31)
If the noise parameter σ is unknown, then we can find the optimal σ∗ value by solving the optimization problem
σ∗ = arg min σ Q (x̂(σ)) (32)
where Q(·) is an image quality penalty such as total variation (TV). In contrast to the SURE, where the σ dependency exists in the loss function in (3) and the neural network must be trained with a new σ value, the optimization problem in (32) is computationally negligible since the precalculated score function estimate l̂′(y) is used regardless of different σ values. A similar on-the-fly parameter estimation can be used for Gamma noise model. As shown in (21), the parameters (α, β) can be adjusted using the precomputed score function estimate l̂′(y). Even for the Poisson noise model, the same trained neural network is used to compute the score function for each input y/ζ . Therefore, our Noise2Score has advantages for the blind application where existing methods cannot be used due to the difficulty of generating samples by changing the noise level parameters without retraining.

heading: 4 Experimental Results
text: Dataset and Implementation detail We evaluated the proposed method for gray-scale and color images in the four benchmark datasets: the gray-scale image dataset contains BSD68 [23] and Set12. The color image dataset contains RGB natural images CBSD68 [23], Kodak dataset. We adopt DIV2K [24], and (C)BSD400 dataset as train data set. We generated the synthetic noise images for each noise distribution. In order to evaluate the proposed method fairly with comparison methods, the same modified U-Net generator [3] is used for all methods. The total epoch was set to 100 and the Adam optimizer [25] was used to train the network. The learning rate was initialized to 2×10−4 for first 50 epoch, and after 50 epoch the learning rate was decayed to 2× 10−5. The proposed method
was implemented in PyTorch [26] with NVidia GeForce GTX 1080-Ti. The network training took about 10 hours. The more detail of implementation are described in Supplementary Material. To deal with blind noise case, as for the quality metric Q(·) we use the TV norm for the Gaussian noises and its combination with data fidelity term for the Poisson and Gamma noise. The details of the quality metric and the discussion on the accuracy of the blind parameter estimation is provided in Supplementary Material.
Gaussian noise Using additive Gaussian noise, our method was compared to (C)BM3D [13], Noise2Void, Noise2Self, SURE, Noise2Noise and supervised learning as shown in Table 2. As expected, the supervised learning and Noise2Noise using multiple noisy realization were the best, but these methods are not practical as explained before. Among the other self-supervised learning approaches, our Noise2Score provide the best performance for all dataset. With the comparison with the supervised learning and Noise2Noise, the performance difference of our method is only margin. Even for blind-noise cases where the noise levels are randomly sampled from Gaussian distribution with σ ∈ [5, 55], Noise2Score with the on-the-fly parameter estimation performed comparable to the supervised learning approach. These results indicate the superior performance of our method. The qualitative comparison in Fig. 2 show that Noise2Score retains the better image details and provides much visually pleasing results compared to other top two self-supervised learning approaches.
Poisson noise For the Poisson noise case, BM3D and SURE were replaced with Anscombe (BM3D+VST) [27] and PURE [12], respectively, as indicated in Table 2, as they are optimized for Possion noise model. Again, our Noise2Score results in significant performance gain compared to other approaches, and the results are comparable to supervised learning approach and Noise2Noise. In addition, to perform experiments when the noise parameters are unknown, we performed experiments with the noise levels randomly sampled between ζ ∈ [0.001,0.1]. We found that our proposed method
produce even higher PSNR than the results of other existing methods even with known parameters, and also provide comparable results to the supervised learning approaches. The qualitative comparison in Fig. 2 with the other two best approach showed that our proposed method provides much cleaner images.
Gamma noise As the extension of BM3D and SURE are not available for Gamma noises, only four comparison methods was adopted for comparison as indicated in Table 2. Here, we set α = β = k. Again, Noise2Score produced best results among self-supervised learning approaches, and provided comparable results to Noise2Noise and supervised learning approaches. Even for the blind noise case, where the noise parameter k randomly selected in k ∈ [40, 120], Noise2Score provided the better PSNR results compared to the other self-supervised learning method even with known parameters, and produces comparable results to supervised learning approaches. The qualitative comparison in Fig. 2 confirm that our Noise2Score provides best reconstruction results.

heading: 5 Conclusion
text: In this work, we provided a novel Bayesian framework for self-supervised image denoising without clean data, which surpasses SURE, PURE, Noise2X, etc. Our novel innovation came from the Tweedie’s formula, which provides explicit representation of denoise images through the score function. By combining with the score-function estimation using AR-DAE, our Noise2Score can be applied to image denoising problem from any exponential family noises. Furthermore, an identical neural network training can be universally used regardless of the noise models, which leads to the noise parameter estimation with minimal complexity. The links to SURE and existing Noise2X were also explained, which clearly showed why our method is a better generalization.
Limitation and negative societal impacts
As a negative societal impact, the failure of image denoising methods could produce side effects. For example, removing both the noise and the texture of the medical images could lead to misdiagnoisis. While the proposed Noise2Score has merits, there are also limitations. In real environments, the prior knowledge of noise distribution may not be available, and the noise model could not be modeled by exponential family noises. Therefore, future extension to such scenario would be beneficial.

heading: Acknowledgment
text: This research was funded by the National Research Foundation (NRF) of Korea grant NRF2020R1A2B5B03001980, ETRI (Electronics and Telecommunications Research Institute)’s internal funds [21YR2500, Development of Digital Biopsy Core Technology for high-precision Diagnosis and Therapy of Senile Disease], and KAIST Key Research Institute (Interdisciplinary Research Group) Project. This work was also supported by Institute of Information & Communications Technology Planning & Evaluation (IITP) grant funded by the Korea government (MSIT) (No.2019-0-00075, Artificial Intelligence Graduate School Program (KAIST)).

Review summaries:
The paper proposes a new approach based on Tweedie's formula to perform self-supervised image denoising. The proposed approach has two steps: (1) learn the score function or gradient of density of the noisy images, and (2) use Tweedie's formula to compute a denoised image using the learned score function. The proposed approach builds on multiple relevant work in the literature, and is well motivated. The model is tested on Gaussian, Poisson and Gamma noise removal tasks on natural images and achieves competitive results. 

This paper considers the problem of image denoising using deep learning. Traditionally, deep denoisers are trained in a supervised way using pairs of noisy and ground-truth images. Recent methods, such as Noise2Noise and those based on SURE, are used to train deep denoisers without ground truth.  This paper proposes Noise2Score as an alternative to existing methods. The idea of Noise2Score is to pre-train a denoising autoencoder (DAE) and use it as an estimator of the score function within the Tweedie’s formula. The key observation is that the pre-training stage of DAE doesn't explicitly need to know the true noise statistics. On other hand, the Tweedie’s formula enables the adaptation of the DAE to specific noise statistics.  The main conceptual contribution of the paper is that it establishes a formal link between the usage of DAEs for approximating score functions and image denoising. The empirical results show the competitive performance of Noise2Score against several well known image denoising methods, including Noise2Noise and SURE.