Planar Nuclear Medicine Images De-Noising Via Wavelet Block Thresholding: a Simulation Study

____________________________________________________________________________________________ ABSTRACT Aims: Wavelet transform is the powerful mathematical tool used for image processing and noise suppression applications in different area of science and engineering. In this technique, selecting optimal threshold for de-noising is still an area of thrust for the researchers. In this paper, we have focused on the de-noising of planar nuclear medicine image using the block thresholding with different block sizes when threshold values vary for Stein, soft and hard thresholdings. Study Design: De-noising of planar nuclear medicine images via wavelet block thresholding Methodology: We simulated planar images of hot region of Carlson phantom by GATE v. 6.1. Noisy image and reference image were produced by imaging time 10 and 40 second, from the hot region of Carlson phantom which is placed next to the simulated gamma camera. Then, we tried to de-noise noisy test image by wavelet transforms and block thresholding methods. For de-noising, we show the evolution of the de-noising the different size of the soft block the results could be obtained by test image which is subjected to the block sizes of 3 and 4. Furthermore, Invariant soft thresholding is found to yield an overly smoothed estimate than orthogonal soft thresholding. Conclusion: Although decreasing of imaging time increases Poisson noise in the acquired nuclear medicine images, using of de-noising technique based on the wavelet transform could improve image degradation, so that the quality of de-noised test image could be compared to the reference image. technique with block size 3 gives the best results for the simulated planar nuclear medicine image in this study.


INTRODUCTION
Nuclear medicine imaging procedures usually have two classifications; the first is those that based on single-photon emitter radionuclides, such as planar and tomographic imaging (single-photon emission computed tomography or SPECT). The second is depending on the two back-to back gamma rays emitted by positron annihilations such as PET (positron emission tomography) imaging.
In planar (or static) single-photon imaging the gamma camera is placed next to the investigated organ, such as thyroid, heart or lungs, and the data are collected for predefined time. Thus, a two-dimensional map of activity distribution is produced.
Although tomographic imaging has rendered planar imaging obsolete for clinical nuclear medicine imaging, but the planar imaging is still important in certain circumstances particularly, for example, for gamma-angiography as well as for bone or renal scans. In the former case, the gamma camera is usually positioned in the 45 degree left anterior oblique projection with 15-30 degrees of caudal tilt. The left anterior oblique angle is then adjusted to acquire clear separation of the ventricles and the caudal tilt can also be adjusted to obtain separation of the atria from the ventricles [1].
Several factors cause poor spatial resolution, low contrast and high noise levels in the acquired nuclear medicine images. These factors could be divided into three parts; first, collimation, uniformity and stability of gamma camera that is related to imaging system, second, factors related to patient, such as positioning and movement, third, parameters related to corrupt recorded events due to different physical interaction of gamma rays and Poisson statistical noise of scintillation camera image where the variance of a measure in a region-of-interest is the sum of the counts [2]. Although increasing of imaging time causes Poisson noise reduction, but this is trouble for the patients and increases patient movements in imaging process. One of the effective ways for enhancement of image quality and noise suppression is use of mathematical procedure such as smoothing filters which strongly affect the quality of the acquired image. These filter functions, which frequently act in Fourier domain, involve a trade-off between the signal-to-noise ratio (SNR) and the spatial resolution of the signal and image processed.
In 1992 a method of signal de-noising depend on adaptive nonlinear filtering in the wavelets transform domain has been presented by D. Donoho [3], because wavelets provide appropriate basis for distinguishing between the image signal and noisy signal. Indeed, multi-resolution representations of signals or images could be produced by the wavelet transforms. This transformation decomposes images into multi scale details. The mother wavelet used in the wavelet transforms are localized around the specified point and has a support size proportional to the definite scale. They are nonzero only over part of the domain represented. Wavelet transform could be preserved sharp transitions in images and signals. This special treatment of edges by wavelet transform is very attractive in image de-noising and processing.
In this paper, we simulated planar images of hot region section of Carlson phantom by GATE v. 6.1 (GEANT4 application for tomographic emission software) [4]. The GATE Monte Carlo package was expanded by the Open-GATE collaboration. It is an open-source extension of the GEANT4 Monte Carlo toolkit and the ROOT object oriented data analysis framework [4]. Noisy image was produced by imaging time 10 second, from hot region section of Carlson phantom which had been placed next to the simulated gamma camera. In general, the only way to reduce the effect of photon noise is to capture more signals. Such that, the ratio of signal to photon noise grows with the square root of the number of photons captured. Thus, reference image was produced by imaging time 40 second as state of art. Then, we tried to de-noise noisy image by wavelet transforms and block thresholding methods. For de-noising, we show the evolution of the de-noising peak signal to noise ratio (PSNR) and Root Mean Square Error (RMSE) when threshold values vary for Stein thresholding, soft thresholding and hard thresholding. And with acquired the best thresholding method with specified parameters, we studied the effect of the block sizes on the PSNR and RMSE of de-noised image. Finally, block thresholding has been applied to a translation invariant wavelet transform and optimal threshold has been calculated for this method. Hereby, several programs were written in MATLAB for de-noising which were based on the MATLAB codes presented by A. Peyré in reverence [5]. In each step, de-noised images were compared with reference image using PSNR and RMSE criterions.

Monte Carlo Simulations
As mentioned, we used GATE v. 6.1 for simulation of planar images. At first, we simulated a gamma camera consists of NaI(Tl) crystal with dimension of 1x19x28 cm 3 . This gamma camera was equipped with a low energy high resolution lead collimator in simulation media. Then, hot regions section of Carlson phantom was simulated and was placed in the next to the simulated gamma camera. This region consist of 9 paired-cylinders as hot regions with 4.7, 5.9, 7.3, 9.2, 11.4, 14.3, 17.9, 22.4 and 38.5 mm diameters in the way of eight pairedholes with same diameter located on the opposite of each other and forming V-shape (as shown in Fig. 1-a). Carlson phantom was frequently filled with water containing 1mCi of Technetium 99-m in practical nuclear medicine study. Then, activities of hot cylinder holes were selected based on specified volume of each cylinder in comparison with the total volume of active water in the Carlson phantom.
At first, the planar image of Carlson phantom was acquired using imaging time 10 second which performed noisy (test) image. For reference image we acquired an image with high counts with increasing of accusation time up to 40 second to acquire a higher quality image for simulated phantom. Fig. 1-b and Fig. 1-c illustrate the image of reference and noisy images, respectively. Then, we transformed these images using a code that was written in MATLAB. This was followed by using noisy initial image along with reference image, where they were normalized and utilized for comparison of results. PSNR and RMSE were utilized for image quality assessment and comparison between produced de-noised images by different methods of wavelet transform.

Wavelet Transform
The general procedure of wavelet de-noising method is that, at first the noisy image is converted into wavelet domain. This transformation cause production of four subbands based on the level of decomposition that frequently shown with "J". Two dimension discrete wavelet transform leads to decomposition of approximate coefficients at level 'J' into four components i.e. the coarse coefficients that constitute the LL subband approximation at level 'J+1' and details subbands in different orientations involving horizontally (HL subband), vertically (LH subband) and diagonally (HH subband). Since the noisy components are dominant in the high frequency spectrum, the probability of existence of noisy components is high in the three higher bands. For each subband (except the low pass or approximation subband) apply proper threshold to the subband coefficients could be smoothed noisy wavelet coefficients. And then, reconstruction the image by employing the inverse discrete wavelet transform may be followed.
In this process, selection of optimal threshold value is important factor for the efficiency of de-noising technique. Threshold value is selected based on different procedure such as Visushrink, Bayes Shrink, Normalshrink, Bivariate shrink and Sure shrink that based on the image and noise priors such as mean and variance [6,7]. Moreover, There are different types of thresholding available consist of Stein, soft and hard thresholding. Fig. 2 shows the one dimensional curves of these thresholding functions.
One of the effective method in the image de-noising via wavelet transform is using of block thresholding technique. Block thresholding takes advantage of the statistical dependency of wavelet coefficients, by computing the attenuation factor on the pre-defined block of wavelet coefficients. When edges and geometric characterize of wavelet transformed images produce collections of high magnitude coefficients, this method especially efficient for denoising [5].
Block de-noising also helps to eliminate artifacts due to large coefficients of isolated noisy in regular areas. The set of wavelet coefficients in each subband is divided into separate blocks, and the block thresholding is applied to the all calculated coefficients within a desired block. The block thresholding strategy can also be applied to wavelet coefficients in translation invariant thresholding. A common method to change a de-noiser into a translation invariant de-noiser is to average the result of translated images [5].

Peak Signal-to-Noise Ratio (PSNR) and Root Mean Square Error (RMSE)
PSNR is an engineering term for the ratio between the maximum possible power of a signal and the power of corrupting noise that affects the fidelity of its representation. PSNR criterion has been used as an estimation of the efficiency of de-noising technique. Where for an original image f and de-noised image f 0 , PSNR calculated as; Where RMSE is the root mean square error between the original (i.e. f 0 ) and the de-noised image (i.e. f) with size M x N: PSNR unit is decibel (dB).

RESULTS AND DISCUSSION
We de-noised the test image by thresholding block of coefficients together. In order to study the function of Stein, soft and hard thresholding, the quality of the de-noised image, which estimated by the PSNR and RMSE criterions, using different threshold values is examined. The noisy planar imaging was used and the threshold value was allowed to span a wide range of numbers in the desired selected block size 4.
We observed that for the given noisy image, the optimal thresholds belong to the soft and Stein thresholding algorithms, as shown in diagrams a and b in Fig. 3, which have maximum PSNR values of 29.49 and 29.26 dB and minimum RMSE values of 1.21 and 1.29, respectively, and the worse results are related to the hard thresholding which has maximum PSNR value of 20.11 dB and minimum RMSE value of 1.63. The PSNR curves related to soft and Stein thresholding, unlike drawn curve by applying hard thresholding, has the same behaviors. This situation also exists for the RMSE curves. Furthermore, the PSNR curve related to soft thresholding reach to its maximum point at the lower value of the threshold value than PSNR curve for Stein thresholding.
Numerically, on the noisy images, soft thresholding gives the best results. Thus, we selected soft thresholding method for de-noising test noisy image which produced the best results in the former PSNR and RMSE diagrams. Then, we applied the different block size, values of 2 to 5, in the block thresholding process for a broad range of threshold values. The obtained results have been shown in Fig. 4. As seen, the best results could be produced by applying block sizes of 3 and 4 in the block thresholding process. Comparing the effect of different size of blocks on the de-noising process of test image by the PSNR and RMSE criterions show that the produced results are very similar to each other, but the best results could be obtained by the test image is subjected to block sizes of 3 and 4, respectively. The best result for block size 3 is about 29.59 dB and for block size 4 is about 29.58 dB which obtains in threshold of 0.028 for both of block sizes. RMSE diagrams for different thresholding method, as shown in Fig. 4-b, show also that the best results could be obtained by threshold value of 0.028 which are about 1.09 and 1.10 for block sizes of 3 and 4, respectively.
Finally, we applied block thresholding to a translation invariant wavelet transform in the denoising process. In order to study the effect of orthogonal and invariant wavelet transform, the quality of the de-noised image, as estimated by the PSNR and RMSE criterions, is examined by use of the various threshold values. The results are shown in Fig. 5. Results show that the translation invariant soft thresholding does a slightly better job than translation orthogonal soft thresholding. So that, the best results for invariant thresholding is about 29.70 dB for PSNR, as shown in Fig.5-a, which is the best result that have been estimated for PSNR of de-noised image, until now. Same situation exists for the RMSE diagram in Fig.  5-b. Thus, this de-noising technique with block size 3 gives the best results for the simulated planar nuclear medicine image in this study.   Fig. 6-b show the images that are subjected to orthogonal and invariant thresholding. The image quality of produced image by applying invariant thresholding seems better than de-noised image by orthogonal thresholding, to some extent. On the other hand, the visual comparison between de-noised image by invariant thresholding, as shown in Fig.  6-a and reference image, as shown in Fig. 6-c, demonstrate the power of block thresholding method with optimal parameters for improvement of quality of noisy test image. So that, the quality of de-noised test image could be compared to the reference image.

CONCLUSION
Diagrams show that soft thresholding has been produced better results than others and then gives more visually pleasant image because of the continuity of soft thresholding function. Hard thresholding is discontinuous and produces abrupt artifacts in the de-noised test image especially when the noise energy is important. Studies are conducted to assess the performance of block thresholding with different block sizes in de-noising of planar nuclear medicine image of the hot region section of Carlson phantom when threshold values vary for Stein, soft and hard thresholding. The results show that soft thresholding with block sizes 3 and 4 removes noise significantly and improves PSNR and RMSE values for produced denoised images. Although. The maximum PSNR and minimum RMSE values related to soft and Stein thresholding are resemble, to some extent.
However, for de-noising images, invariant soft thresholding is found to yield an overly smoothed estimate than orthogonal soft thresholding as seen in Fig.s 5 and 6. Indeed, the PSNR improvement could be produced by the translation invariant thresholding, because averaging process in the translation invariant de-noiser technique causes reduction of oscillating artifacts of orthogonal thresholding.