IOP PUBLISHING
PHYSICS IN MEDICINE AND BIOLOGY
Phys. Med. Biol. 56 (2011) 5221–5234
doi:10.1088/0031-9155/56/16/009
Maximum likelihood estimation-based denoising of magnetic resonance images using restricted local neighborhoods Jeny Rajan1 , Ben Jeurissen1 , Marleen Verhoye2 , Johan Van Audekerke2 and Jan Sijbers1 1 2
Vision Lab, University of Antwerp, Belgium Bio-Imaging Lab, University of Antwerp, Belgium
E-mail:
[email protected]
Received 31 March 2011, in final form 23 June 2011 Published 26 July 2011 Online at stacks.iop.org/PMB/56/5221 Abstract In this paper, we propose a method to denoise magnitude magnetic resonance (MR) images, which are Rician distributed. Conventionally, maximum likelihood methods incorporate the Rice distribution to estimate the true, underlying signal from a local neighborhood within which the signal is assumed to be constant. However, if this assumption is not met, such filtering will lead to blurred edges and loss of fine structures. As a solution to this problem, we put forward the concept of restricted local neighborhoods where the true intensity for each noisy pixel is estimated from a set of preselected neighboring pixels. To this end, a reference image is created from the noisy image using a recently proposed nonlocal means algorithm. This reference image is used as a prior for further noise reduction. A scheme is developed to locally select an appropriate subset of pixels from which the underlying signal is estimated. Experimental results based on the peak signal to noise ratio, structural similarity index matrix, Bhattacharyya coefficient and mean absolute difference from synthetic and real MR images demonstrate the superior performance of the proposed method over other state-of-the-art methods. (Some figures in this article are in colour only in the electronic version)
1. Introduction Data acquired by an MRI system are inherently corrupted by statistical noise of which the primary sources are electronic and dielectric and inductive coupling to the conducting solution inside the body (C´ardens-Blanco et al 2008). Noise remains one of the main causes of quality deterioration in MRI and is a subject in a large number of papers in MRI literature, e.g. 0031-9155/11/165221+14$33.00
© 2011 Institute of Physics and Engineering in Medicine
Printed in the UK
5221
5222
J Rajan et al
Edelstein et al (1984), Henkelman (1985), McVeigh et al (1985), Bernstein et al (1989), Gudbjartsson and Patz (1995), Macovski (1996). Other than visual analysis, processing techniques such as segmentation, registration or tensor estimation in diffusion tensor MRI (DT-MRI) will be affected or biased due to noise (Jones and Basser 2004, Aja-Fern´andez et al 2010). Noise can be naturally minimized by averaging images after multiple acquisitions. This, however, may not be feasible in clinical and small animal MRI where there is an increasing need for speed (Manj´on et al 2008). Thus, post-processing techniques to remove noise in the acquired data are important. Several filtering techniques to improve the quality of the magnetic resonance (MR) images have been proposed in the literature. Most of the methods proposed earlier can be mainly classified as either based on partial differential equations (PDEs), wavelets or nonlocal means (NLM). A PDE-based approach for filtering MRI was first attempted by Gerig et al (1992). In their work, they demonstrated that anisotropic diffusion is an effective filtering technique for MRI in the sense that it can significantly decrease the image noise and simultaneously preserve fine details in the image. A major drawback of their method, however, was the incorrect assumption about the noise distribution. The noise was assumed to be Gaussian instead of Rician, as a result of which a bias is introduced in the filtered image. Such a bias becomes particularly important in low SNR MR images, such as diffusion weighted images (AjaFern´andez et al 2008). To account for the Rice distribution, an adaptive anisotropic diffusion method for magnitude MR data was proposed by Sijbers et al (1999). Finally, Samsonov and Johnson (2004) presented a noise adaptive nonlinear diffusion filtering technique to denoise MR images with spatially varying noise levels. All aforementioned PDE methods are based on classical second order Perona–Malik (Perona and Malik 1990) anisotropic diffusion. Although such methods are effective in denoising images, they tend to cause staircase effects in the filtered images (You and Kaveh 2000). To reduce this effect, a noise removal algorithm for MRI based on fourth order PDE was suggested by Lysaker et al (2003). The main strength of this method is its ability to process signals with a smooth change in the intensity value. Basu et al (2006) used a data likelihood term combined with Perona–Malik anisotropic diffusion to effectively denoise an MR image. Recently Krissian and Aja-Fern´andez (2009) proposed a noise-driven anisotropic diffusion filter for denoising MR images, in which the diffusion is controlled by the local statistics in the image derived from the linear minimum mean square error (LMMSE) estimator for the Rician model. A second class of noise filtering schemes are wavelet based (Nowak 1999, Wood and Johnson 1999, Pizurica et al 2003). These algorithms exploit the decorrelating properties of the wavelet transform to suppress noise coefficients using statistical inference. Among these methods, Nowak (1999) got much attention, in which a bias removal was proposed, based on the observation that the squared Rician distributed data exhibits a chi-squared distribution with two degrees of freedom. However, in Manj´on et al (2008), it is mentioned that the aforementioned wavelet-based filters may introduce characteristic artifacts that can be quite problematic. A trilateral filter was proposed in Wong and Chung (2004) to take into account the local structure in the image, in addition to intensity and geometric features. Recently, Delakis et al (2007) proposed a wavelet-based denoising algorithm for images acquired with parallel MRI. During the past few years, NLM-based denoising methods have gained much popularity (Buades et al 2005). Manj´on et al (2008) were the first to attempt MRI denoising with the NLM approach. Coup´e et al (2008) proposed an optimized blockwise version of the NLM algorithm for denoising MR images. Wiest-Daessl´e et al (2008) suggested an adaptive NLM
Maximum likelihood estimation-based denoising of magnetic resonance images
5223
filter which behave better for Rician corrupted MR data. Aja-Fern´andez et al (Aja-Fern´andez and Krissian 2008) incorporated the LMMSE approach in NLM to make it adaptive for Rician data. Recently, Manj´on et al (2010) proposed an adaptive NLM denoising method for MR images with spatially varying noise levels. Apart from the above-discussed methods, other popular approaches proposed for MR denoising are the methods based on the maximum likelihood (ML) estimation (Sijbers and den Dekker 2004, He and Greenshields 2009) and the filtering based on smoothing splines (Koay et al 2009). In this paper, we propose a method to denoise magnitude MRI based on the ML estimation method using a restricted local neighborhood. ML estimation can be applied locally (also referred as local ML (LML)) or nonlocally as proposed by He and Greenshields (2009), in which ML estimation is applied on a set of pixels selected based on the similarity of the neighborhood. One drawback of the nonlocal ML (NLML) estimation method is the use of a fixed sample size for the ML estimation, which can cause either under- or over-smoothing. The disadvantage of the LML method is the blurring of edges and the distortion of fine structures in the image (He and Greenshields 2009). This is because the assumption that the signal in the selected small neighborhood is constant is generally not valid. In this paper, we put forward the concept of a restricted local neighborhood as a solution to this problem. This paper is organized as follows. Section 2 gives a short overview about the noise characteristics in MRI. Section 3 elaborates the proposed method. Section 4 presents the experimental results, comparative evaluation and discussion. Finally, conclusions and remarks are drawn in section 5. 2. Noise characteristics in MRI The raw, complex MR data acquired in the Fourier domain are characterized by a zero mean Gaussian probability density function (PDF). After the inverse Fourier transform, the noise distribution in the real and imaginary components will still be Gaussian due to the linearity and the orthogonality of the Fourier transform. However, due to the subsequent transform to a magnitude image, the data will no longer be Gaussian but Rician distributed. For an MR magnitude image defined on a discrete grid , M = {mi |i ∈ }, the probability distribution function of mi, with underlying signal A, is (Rice 1945, Bernstein et al 1989) 2 m2 i +A mi − 2σ Am i 2 g I p(mi |A, σg ) = 2 e (1) (mi ), 0 σg σg2 where I0 (.) is the 0th order modified Bessel function of the first kind, (.) is the Heaviside step function and σg2 denotes the variance of the Gaussian noise in the complex MR data. If A equals zero, the Rician PDF simplifies to a Rayleigh distribution (Edelstein et al 1984) m2
p(mi |σg ) =
mi − 2σig2 e (mi ). σg2
(2)
At high SNR, the Rician PDF approaches a Gaussian PDF with mean A and variance σg2 : p(mi |σg ) =
1 2π σg2
−
e
(mi −A)2 2σg2
(mi ),
(3)
i.e. Rician noise in magnitude MR images assumes a Gaussian distribution when the SNR is high and Rayleigh distribution for low SNR. Proper estimation of the noise variance, σg2 , is important for effective denoising of MRI. Many methods have been proposed in the literature for the estimation of the noise level
5224
J Rajan et al
from MRI. A survey of these methods is given in Aja-Fern´andez et al (2009). Most of the noise estimation methods for MRI estimate the noise level from the background area of the magnitude MR image, which is known to be Rayleigh distributed (Sijbers et al 2007, AjaFern´andez et al 2008). These methods need a certain amount of background pixels to perform proper estimation and can fail in the absence of the background region (Aja-Fern´andez et al 2010). Recently, in Rajan et al (2010), Coup´e et al (2010) and Manj´on et al (2010), methods are proposed to estimate noise in the absence of background. In Rajan et al (2010), two different approaches based on the local estimation of variance using ML and skewness are proposed and the methods suggested in Coup´e et al (2010) and Manj´on et al (2010) are based on the correction factor suggested in Koay and Basser (2006). In this work, we followed the methods in Rajan et al (2010) for noise variance estimation. 3. Methods 3.1. Signal estimation using LML Let m1 , m2 , . . . , mn be n statistically independent observations from a region of constant signal intensity A. Then, the joint pdf of the observations is 2 n m2 i +A mi − 2σ Am i 2 g I e p({mi }|A, σg ) = . (4) 0 σ2 σg2 i=1 g The ML estimate of A can now be computed by maximizing the likelihood function L(A) or equivalently ln L(A), with respect to A (Sijbers and den Dekker 2004): n n n m2i + A2 mi Ami + ln I0 ln L = ln − (5) σg2 2σg2 σg2 i=1 i=1 i=1 and AML = arg{max(ln L)}. A
(6)
In a region of constant signal amplitude, the true underlying signal can be estimated from the noisy signal using (6). The straightforward approach to denoise MR images using the above-explained ML estimation method is to apply the method locally for each pixel with the assumption that the underlying area is constant for a small region. However, this assumption is not valid for regions with edges and fine structures and as a result these regions will get blurred. This effect can be observed from figure 1. It can also be observed from figure 1(c) and figure 1(d) that the image denoised with a window size of 3 × 3 × 3 is, not surprisingly, less blurred than the one denoised with a window size of 5 × 5 × 3. However, selecting a very small window size is less effective for denoising, especially in smooth areas, since the number of samples for estimation is smaller. In the proposed method, we consider only the pixels that have an underlying gray-level value close to that of the center pixel in the local window for the true signal estimation, instead of selecting all pixels in the window. This approach can reduce the side effects of direct local ML estimation. However, selection of pixels with similar underlying gray value from a noisy image is a difficult problem. To solve this issue, we used a reference image. The reference image is created by applying the NLM algorithm over the noisy image. Based on the information from the reference image, the corresponding pixels are selected from the noisy image for true signal estimation.
Maximum likelihood estimation-based denoising of magnetic resonance images
(a)
(b)
(c)
5225
(d)
Figure 1. Denoising MRI with the LML method: (a) original image and (b) original image corrupted with Rician noise of σg = 15. (c) Noisy image denoised with LML of window size 3 × 3 × 3. (d) Noisy image denoised with LML of window size 5 × 5 × 3. The blurring effect generated by the LML method can be observed from the denoised images.
3.2. Nonlocal means algorithm The NLM method which we used to create the reference image is briefly explained in this section. The NLM method was proposed by Buades et al (2005) and is based on the Markovian hypothesis which states that pixels with a similar neighborhood have a similar gray level value. Given a noisy image M = {mi |i ∈ }, where represents the image space and mi corresponds to the noisy image value at location i, the filtered value at a point i, υi , is calculated using the NLM method as a weighted average of all the pixels in the image (Buades et al 2005): 1 wi,j mj , (7) υi = Ci j s
where s represents the neighborhood pixels and Ci = wi,j
(8)
j s
is a normalization constant and the weight wi,j is determined by the similarity of the Gaussian neighborhood between pixels i and j , which can be expressed as Ni − Nj 22,a wi,j = exp − , (9) h2
where Ni denotes a square neighborhood centered at pixel i, · 22,a is a Gaussian-weighted Euclidean distance function, a is the standard deviation of the Gaussian kernel and h acts as a degree of filtering. 3.3. Signal estimation using restricted LML To overcome the drawback of the local ML (LML) estimation method, we propose a restricted LML (RLML) estimation method. In RLML, only the pixels in the local neighborhood of the noisy pixel mi that have an underlying gray level value close to the underlying gray level value of mi, will be considered for the true signal estimation. However, as mentioned earlier, selection of pixels with similar underlying gray value from a noisy image is a difficult problem. To this end, we create a reference image using the above-mentioned NLM method. Now, to denoise a noisy pixel mi at i, a list li is created from the neighbors of mi: li = {mj , (j ∈ m )| abs(f (mj ) − f (mi )) < t},
(10)
5226
J Rajan et al 100
3 x 3 neighborhood
1500
80
5 x 5 neighborhood
60
1000
40 500 20 0 100 110 120 130 140 150 160 170 180 190 200
0
0
20
40
60
80 100 120 140 160 180 200 220 240
intensity level
distribution of local range
(b)
(c)
(a)
Figure 2. Computation of threshold from ROI and from the mode of the local distribution of range: (a) reference image with ROI (marked as red), (b) distribution of intensity in the ROI, (c) distribution of local range computed using 3 × 3 and 5 × 5 neighborhood.
where m represents the neighborhood space around mj, f (mj ) = υj and f (mi ) = υi . The threshold t used for the classification is calculated from the reference image as the range of the intensity dispersion of a uniform area. This can also be automatically computed from the reference image, υ, as the mode of all the local distributions of the range computed around the neighborhood of each pixel: t = mode{range(υF )w }
(11)
where υF represents the foreground region of the reference image and w is the neighborhood window size. Contrary to complex-valued images, where the noise is Gaussian distributed, in magnitude images the range of intensity dispersion will depend on the local SNR. Hence, to reduce the error in the classification of pixels in the reference image, the threshold t is computed only from the foreground region of the image. Figure 2 shows the distribution of the local range for a uniform region in the reference image and also the distribution of the local range computed with a neighborhood size of 3 × 3 and 5 × 5. It can be observed from the plots in the figure that the mode of the distributions of the local range is close to the actual range of the intensity dispersion of the selected homogeneous region. Once a list is created as mentioned in (10), the denoised pixel Ai at location i can be computed by substituting the values in the list li as the Rician-distributed magnitude data points in (5) and then maximizing the log likelihood function. Applying this procedure to all pixels in the noisy image will give the denoised image. The method is summarized in algorithm 1. 3.4. Spatially varying noise levels When measuring the noise in MR images, it is often assumed that the noise is not spatially varying across the image. Parallel MRI (pMRI) acquisition techniques such as sensitivity encoding (SENSE) or generalized autocalibrating partially parallel acquisitions (GRAPPA) introduce a spatially varying noise variance across the image (Manj´on et al 2010, AjaFern´andez et al 2010). To deal with spatially varying noise, local noise variance estimation should be introduced (Manj´on et al 2010). Therefore, we estimate both the noise variance and the underlying true intensity simultaneously by maximizing the log likelihood function in (5) with respect to both A and σg2 . The ML estimate is then found from the global maximum of ln L w.r.t. A and σg2 (Sijbers and den Dekker 2004): AML , σg2 ML = arg{max(ln L)}. (12) A,σg2
Maximum likelihood estimation-based denoising of magnetic resonance images
5227
Algorithm 1. Algorithm for signal estimation using RLML 1: Estimate the noise standard deviation σg from the input magnitude image M using the method described in Rajan et al (2010) 2: Create the reference image υ using (7) 3: Compute the threshold t from υ by applying (11) using a neighborhood window of size n × n × n 4: for every pixel mi of M do 5: Create the list li as mentioned in (10) using a neighborhood window of size w × w × w 6: Substitute the values in the list li in (5) 7: Estimate Ai by maximizing (5) with respect to the unknown true intensity 8: end for
4. Experiments and results To evaluate and compare the proposed methods with state-of-the-art methods, we did experiments on both synthetic and real MR images. To conduct the experiments on synthetic data, we used the standard MR image phantom of the brain obtained from the Brainweb database (Cocosco et al 1997) (T1-weighted, intensity values in the range 0–255) and for the experiments on DWI, we simulated a set of DWI of the human brain using the methods in Leemans et al (2005, 2009). For the experiments on real data, we used the MR images of a kiwi fruit. The proposed algorithm was compared with the following recently proposed methods. i. NLML: the nonlocal ML method (He and Greenshields 2009) with search window size = 11 × 11 × 11, neighborhood size = 3 × 3 × 3 and sample size 25. ii. UNLM: unbaised nonlocal means (Manj´on et al 2008) with search window size = 11 × 11 × 11, neighborhood size = 3 × 3 × 3 and the value of the decay parameter h = σg . iii. RNRAD: noise-driven anisotropic diffusion filter for MRI (Krissian and Aja-Fern´andez 2009). Local statistics were computed on a 3 × 3 × 3 neighborhood and time step as dt = 1/6. For simulations, the number of iterations was chosen as the one with best PSNR when compared with the original image. iv. ARNLM: adaptive Rician nonlocal means filter with wavelet mixing (Manj´on et al 2010). For a quantitative analysis of the denoising methods, we used the peak signal to noise ratio (PSNR), the structural similarity index matrix (SSIM) (Wang et al 2004), Bhattacharyya coefficient (BC) (Bhattacharyya 1943) and the mean absolute difference (MAD). Figure 3 shows the visual quality comparison of the image denoised with NLML, UNLM, RNRAD and the proposed RLML method. This experiment was conducted on a 2D slice of the synthetic image of the brain in the 3D environment after corrupting the image by noise with σg = 20. The proposed RLML filter was executed with the following parameters: neighborhood size for denoising as 7 × 7 × 3 and the neighborhood size for the local computation of the range as 3 × 3 × 3. In visual analysis, the expectations are as follows: (i) perceptually flat regions should be smooth as possible, (ii) image edges and corners should be well preserved, (iii) texture detail should not be lost and (iv) few or ideally no artifacts (He and Greenshields 2009, Chen et al 2010). It can be observed from figure 3 that the image denoised with the proposed method is
5228
J Rajan et al
(a)
(b)
(c)
(d)
(e)
(f )
Figure 3. Denoising MRI with various methods: (a) original image, (b) original image corrupted with Rician noise of σg = 20, (c) part (b) denoised with the NLML method, (d) part (b) denoised with the UNLM method, (e) part (b) denoised with the RNRAD method and (f) part (b) denoised with RLML method.
closer to the original one (based on the above-mentioned criteria) than the images denoised with the other approaches. Figure 4 shows the quantitative analysis of the proposed method with other recently proposed methods. All the methods to which the proposed method was compared were based on the Rician noise model. In the quantitative analysis, the background was excluded, that is, only the area of the image inside the skull was considered. It can be studied from the plots in figure 4 (taken as a mean of 15 experiments) that the proposed method outperforms other approaches in terms of PSNR, SSIM, BC and MAD. Figures 5 and 6 show the results of the experiment conducted on the synthetic image of the brain with spatially varying noise which was generated using the method mentioned in Manj´on et al (2010). The noisy images are then denoised with the proposed RLML and also with the recently proposed ARNLM filter. Since the NLM approach is not effective in the case of spatially varying noise, we used the adaptive NLM (ANLM) proposed in Manj´on et al (2010) as the reference image for threshold computation. Both visual and quantitative analyses show that the image processed with RLML is more effective than ARNLM in denoising spatially varying noise. For the experiments on DW-MRI, we simulated a set of DW images of the human brain with the following parameters: b = 1200 s mm−2 , voxel size = 1.79 × 1.79 × 2.4 mm3 ,
Maximum likelihood estimation-based denoising of magnetic resonance images
5229
1
40
Noisy NLML UNLM RNRAD RLML
35
0.9 0.8
Mean SSIM
PSNR
30
25
0.7 0.6 0.5
20
Noisy NLML UNLM RNRAD RLML
0.4 15 0.3 10
0
10
20
30
40
50
60
70
0.2
80
0
10
20
30
σ
(a)
60
70
80
50
60
70
80
(b) 70
0.95
Noisy NLML UNLM RNRAD RLML
60
0.9 50 0.85 0.8
MAD
BC
50
g
1
0.75
40
30
0.7
Noisy NLML UNLM RNRAD RLML
0.65 0.6 0.55
40
σ
g
0
10
20
20
10
30
40
50
60
70
80
0
0
10
20
30
40
σ
σ
(c)
(d)
g
g
Figure 4. Quantitative analysis of the proposed method with other recently proposed methods based on (a) PSNR, (b) SSIM, (c) BC and (d) MAD for image corrupted with Rician noise of σg varying from 5 to 80.
(a)
(b)
(c)
(d)
Figure 5. Denoising MRI with spatially varying noise levels: (a) noisy image with spatially varying noise, (b) noise modulation map, (c) noisy image denoised with the ARNLM method and (d) noisy image denoised with the proposed RLML method.
image size = 107 × 79 × 60 and the gradient orientations = 15. The DW images were then corrupted with Rician noise of σg = 100. The denoised DW images were then generated using the proposed method and other methods mentioned earlier. The fractional anisotropy (FA) map was calculated for the ground truth, the noisy and the denoised methods
5230
J Rajan et al 1
40
Noisy ARNLM RLML
35
0.8
Mean SSIM
PSNR
30
25
20
15
0.7 0.6 0.5 0.4
10
5
Noisy ARNLM RLML
0.9
0.3
0
10
20
30
40
50
60
70
0.2
80
0
10
20
30
σ
40
50
60
70
80
50
60
70
80
σ
g
g
(a)
(b) 120
0.95
Noisy ARNLM RLML
0.9 100
0.85 0.8
80
MAD
BC
0.75 0.7
60
0.65 40
0.6 0.55
Noisy ARNLM RLML
0.5 0.45
0
10
20
20
30
40
50
60
70
80
0
0
10
20
30
40
σg
σg
(c)
(d)
Figure 6. Performance comparison of the RLML method against a recently proposed ARNLM method in terms of (a) PSNR, (b) SSIM, (c) BC and (d) MAD for image corrupted with spatially varying Rician noise of σg varying from 5 to 80.
(see figure 7). Figure 8 shows the absolute difference of the estimated FA map with the original FA map for the noisy and all denoising methods. The MAD of the FA residuals shows that the error in the FA map computed from the image denoised with the proposed method is comparatively less than the other methods. Figure 9 shows the results of applying the different denoising methods on an MR image of a kiwi fruit. Two sets of kiwi fruit images were reconstructed, one without averaging and the other by averaging 12 acquisitions. Averaging was done in the complex k-space. The denoising algorithms were then applied over the image reconstructed without averaging and the resultant denoised image was compared with the image reconstructed by averaging 12 acquisitions. It can be observed from the images that the visual results are much better for RLML in terms of image contrast. Quantitative analysis of the experiments on the kiwi fruit was done based on the second order moment of the Rice distribution: E[M 2 ] = A2 + 2σg2 .
(13)
Maximum likelihood estimation-based denoising of magnetic resonance images
5231
(a)
(b)
(c)
(d)
(e)
(f )
Figure 7. Experiments on the DWI atlas of the human brain: (a) FA map computed from the ground truth (b) FA map computed after corrupting the original image with Rician noise of σg = 100, (c) FA map computed from the image denoised with the NLML method, (d) FA map computed from the image denoised with the UNLM method, (e) FA map computed from the image denoised with the RNRAD method, (f) FA map computed from the image denoised with the proposed RLML method.
(a)
(b)
(c)
(d)
(e)
Figure 8. Absolute FA residuals: (a) noisy (MAD: 0.1189), (b) NLML method (MAD: 0.0560), (c) UNLM method (MAD: 0.0531), (d) RNRAD method (MAD: 0.0640), (e) proposed RLML method (MAD: 0.0505).
5232
J Rajan et al
(a)
(b)
(c)
(d)
(e)
(f )
Figure 9. Experiments on the MR image of a kiwi fruit: (a) original image reconstructed with 1 average, (b) original image reconstructed with 12 averages, (c) part (a) denoised with the NLML method (sample size :18), (d) part (a) denoised with the UNLM method, (e) part (a) denoised with the RNRAD method and (f) part (a) denoised with the RLML method. Table 1. Quantitative analysis of the proposed method with other recently proposed methods based on (14). This experiment was conducted on the MR image of a kiwi fruit.
Method Kiwi fruit with 1 average (noise standard deviation = 27.5) Kiwi fruit with 12 averages NLML UNLM RNRAD RLML
Estimated σg based on (14) 31.01 32.84 38.17 23.75 30.18
If we assume the denoised image A as the ground truth, then σg2 can be estimated as
M 2 − A2 , (14) σg2 = 2 where denotes the spatial average of the whole image. The closer the estimated σg2 to the actual noise variance, the closer the denoised image A to the ground truth. The result of this experiment is given in table 1. It is clear from the table that the σg estimated for the image denoised with the proposed RLML method is more close to the actual noise variance. The
Maximum likelihood estimation-based denoising of magnetic resonance images
5233
estimated standard deviation of the noise from the image acquired without averaging is 27.5. The expected value of σg (estimated using (14)) after reconstructing an image by averaging multiple acquisitions or by denoising is a value close to 27.5. Among the denoising methods considered, the estimated value of σg more close to 27.5 is the one denoised with RLML. This experiment on the real data set additionally indicates that the image denoised with the proposed method is more close to the ground truth than other methods. 5. Conclusion A new method to denoise the MR images by applying the ML method locally to restricted neighborhood is proposed in this paper. A scheme is developed to locally select the appropriate subset of pixels from the neighborhood of each pixel. Through this approach, the side effects of the LML method, the blurring effect and the distortion of fine structures can be reduced. Experiments have been carried out on simulated and real data sets. Quantitative analysis at various noise levels based on the similarity measures, PSNR, SSIM, BC and MAD shows that the proposed method is more effective than other state-of-the-art methods. Experiments were also performed on DW images to prove the efficacy of the proposed method. The MAD of the FA residuals shows that the image denoised with the proposed method is more close to the ground truth. Acknowledgments This work was financially supported by the Inter-University Attraction Poles Program 6-38 of the Belgian Science Policy, and by the SBO-project QUANTIVIAM (060819) of the Institute for the Promotion of Innovation through Science and Technology in Flanders (IWT-Vlaanderen). References Aja-Fern´andez S et al 2008 Noise and signal estimation in magnitude MRI and Rician distributed images: a LMMSE approach IEEE Trans. Image Process. 17 1383–98 Aja-Fern´andez S et al 2009 Noise estimation in single- and multiple-coil magnetic resonance data based on statistical models Magn. Reson. Imaging 27 1397–409 Aja-Fern´andez S et al 2010 About the background distribution in MR data: a local variance study Magn. Reson. Imaging 28 739–52 Aja-Fern´andez S and Krissian K 2008 An unbiased non-local means scheme for DWI filtering Proc. MICCAI Workshop on Computational Diffusion MRI pp 277–84 Basu S et al 2006 Rician noise removal in diffusion tensor MRI Proc. MICCAI pp 117–25 Bernstein M A et al 1989 Improved detectability in low signal-to-noise ratio magnetic resonance images by means of phase-corrected real construction Med. Phys. 16 813–7 Bhattacharyya A 1943 On the measure of divergence between two statistical populations defined by their probability distributions Bull. Calcutta Math. Soc. 35 99–109 Buades A et al 2005 A review of image denoising algorithms, with a new one Multiscale Model. Simul. 4 490–530 C´ardens-Blanco A et al 2008 Noise in magnitude magnetic resonance images Concepts Magn. Reson. A 32 409–16 Chen Q et al 2010 Homogenity similarity based image denoising Pattern Recognit. 43 4089–100 Cocosco C A et al 1997 Brainweb: online interface to a 3D MRI simulated brain database NeuroImage 5 S425 http://www.bic.mni.mcgill.ca/brainweb/ Coup´e P et al 2008 An optimized blockwise nonlocal means denoising filter for 3D magnetic resonance images IEEE Trans. Med. Imaging 27 425–41 Coup´e P et al 2010 Robust Rician noise estimation for MR images Med. Image Anal. 14 483–93 Delakis I et al 2007 Wavelet based denoising algorithm for images acquired with parallel magnetic resonance imaging (MRI) Phys. Med. Biol. 52 3741–51
5234
J Rajan et al
Edelstein W A et al 1984 A signal-to-noise calibration procedure for NMR imaging systems Med. Phys. 11 180–5 Gerig G et al 1992 Nonlinear anisotropic filtering of MRI data IEEE Trans. Med. Imaging 11 221–32 Gudbjartsson H and Patz S 1995 The Rician distribution of noisy MRI data Magn. Reson. Med. 34 910–4 He L and Greenshields I R 2009 A nonlocal maximum likelihood estimation method for Rician noise reduction in MR images IEEE Trans. Med. Imaging 28 165–72 Henkelman R M 1985 Measurement of signal intensities in the presence of noise in MR images Med. Phys. 12 232–3 Jones D and Basser P 2004 ‘Squashing peanuts and smashing pumpkins’: how noise distorts diffusion MR data Magn. Reson. Med. 52 979–93 Koay C and Basser P 2006 Analytically exact correction scheme for signal extraction from noisy magnitude MR signals J. Magn. Reson. 179 317–22 Koay C et al 2009 A signal transformational framework for breaking the noise floor and its applications in MRI J. Magn. Reson. 197 108–19 Krissian K and Aja-Fern´andez S 2009 Noise driven anisotropic diffusion filtering of MRI IEEE Trans. Image Process. 18 2265–74 Leemans A et al 2005 Mathematical framework for simulating diffusion tensor MR neural fiber bundles Magn. Reson. Med. 53 944–53 Leemans A et al 2009 ExploreDTI: a graphical toolbox for processing, analyzing, and visualizing diffusion MR data Proc. Intl Soc. Mag. Reson. Med. 3537 Lysaker M et al 2003 Noise removal using fourth-order partial differential equation with applications to medical magnetic resonance images in space and time IEEE Trans. Image Process. 12 1579–90 Macovski A 1996 Noise in MRI Magn. Reson. Med. 36 494–7 Manj´on J V et al 2008 MRI denoising using non local means Med. Image Anal. 12 514–23 Manj´on J V et al 2010 Adaptive non local means denoising of MR images with spatially varying noise levels J. Magn. Reson. Imaging 31 192–203 McVeigh E R et al 1985 Noise and filtration in magnetic resonance imaging Med. Phys. 12 586–91 Nowak R D 1999 Wavelet based Rician noise removal for magnetic resonance images IEEE Trans. Image Process. 10 1408–19 Perona P and Malik J 1990 Scale space and edge detection using anisotropic diffusion IEEE Trans. Pattern Anal. Machine Intell. 12 629–39 Pizurica A et al 2003 A versatile wavelet domain filteration technique for medical imaging IEEE Trans. Med. Imaging 22 323–31 Rajan J et al 2010 Noise measurement from magnitude MRI using local estimates of variance and skewness Phys. Med. Biol. 55 N441–9 Rice S O 1945 Mathematical analysis of random noise Bell Syst. Tech. J. 24 46–156 Samsonov A A and Johnson C 2004 Noise adaptive nonlinear diffusion filtering of MR images with spatially varying noise levels Magn. Reson. Med. 52 798–806 Sijbers J and den Dekker A J 2004 Maximum likelihood estimation of signal amplitude and noise variance from MR data Magn. Reson. Med. 51 586–94 Sijbers J et al 1999 Adaptive anisotropic noise filtering for magnitude MR data Magn. Reson. Imaging 17 1533–9 Sijbers J et al 2007 Automatic estimation of the noise variance from the histogram of a magnetic resonance image Phys. Med. Biol. 52 1335–48 Wang Z et al 2004 Image quality assessment: from error visibility to structural similarity IEEE Trans. Image Process. 13 600–12 Wiest-Daessl´e N et al 2008 Rician noise removal by non-local means filtering for low signal-to-noise ratio MRI: applications to DT-MRI Proc. MICCAI pp 171–9 Wong W and Chung A 2004 Trilateral filtering: a nonlinear noise reduction technique for MRI Proc. Intl. Soc. Mag. Reson. Med. 2218 Wood J C and Johnson K M 1999 Wavelet packet denoising of magnetic resonance images: importance of Rician noise at low SNR Magn. Reson. Med. 41 631–5 You Y and Kaveh M 2000 Fourth order partial differential equations for noise removal IEEE Trans. Image Process. 9 1723–30
