A Log-Exponential Shrinkage Technique for DeNoising of Ultrasound Images for Breast Cancer

Ultrasound in diagnostic imaging is well known for its safety and accessibility. But its efficiency for diagnosis is always limited by the presence of noise. So, in this study, a Log-Exponential shrinkage technique is presented for denoising of ultrasound images. A Combinational filter was designed for the removal of additive noise without losing any details. The speckle noise after homomorphic transformation follows Gaussian distribution and the conventional median estimator has very low accuracy for Gaussian distribution. The scale parameter calculated from the sub-band coefficients after homomorphic transformation was utilized to design the estimator. For shrinkage of wavelet coefficients, a multi-scale thresholding function was designed, with better flexibility. The proposed technique was tested for both medical and standard images. A significant improvement was observed in the estimation of speckle noise variance. For quantitative evaluation of the proposed technique with existing denoising methods, Mean Squared Error (MSE), Structural Similarity Index (SSIM), and Peak Signal to Noise Ratio (PSNR) were used. At the highest noise variance, the minimum improvement achieved by the proposed denoising technique in PSNR, SSIM, and MSE was 10.65%, 23.21%, and 30.46% respectively.


Introduction
Breast cancer has resulted in an immense count of casualties among women [1] and numbers are still rising. Developing countries like India is now observing an increase in the count of young female patients [2]. Early detection has an admirable influence on the mortality rate. Out of many imaging techniques in practice for cancer diagnosis, one involving Ultrasound is always preferred. Along with its proficiency in non-invasive diagnosis, it is very economic and causes zero tissue damage because of no harmful radiation involved, and so allows frequent diagnosis [3], [4]. Dispersed scattering of ultrasound waves results in non-coherent echoes. A locally correlated and granular texture noise appears because of these out-of-phase echoes [5]. They have both multiplicative (called Speckle noise) as well as additive effects [6].
The Major limitation in the application of spatial filters is the choice of kernel size, a small kernel leads to less efficient denoising while blurring in image increases with the increasing kernel size. In wavelet based denoising methodology the noise variance estimation is a preeminent factor for threshold calculation. Initially, D. L. Donoho and I. M. Johnstone [21] applied the Median Absolute Deviation (MAD) based estimator [24] (called as MAD in rest of the study)for this, afterward nearly every study adopted the same estimator. The estimator is recommended for completely symmetrical distribution but its estimation differs a lot from the original value in the case of Gaussian distribution [25]. The speckle noise follows a Rayleigh distribution of finite variance [26], but after homomorphic transformation it follows Gaussian distribution [27]. Moreover, the soft thresholding diminishes a significant number of pixels and leads to undesired dark spots reducing the image quality for further processing.
In this study, an adaptive multi-scale wavelet thresholding with a novel noise variance estimator for denoising of ultrasound images is proposed. To remove the zero mean additive noise, smoothing technique like the Gaussian filter was applied. To avoid the blurring of edges and preserving the details, the gradient of the input image was also used. A novel technique for the estimation of speckle noise variance has been introduced. The proposed estimator is independent of the symmetrical assumption of the distribution. An adaptive thresholding operation with a unique function was also proposed in this study, that provides more flexibility to thresholding operation.

Methodology Datasets
In this study for breast cancer ultrasound images, two different datasets have been utilized. Both of them are published online and are available for researchers free of cost.  [29].
To further authenticate the applicability of proposed technique, liver ultrasound [30], Lena and Barbara images [31] were also utilized in this study. Note that all datasets are published after being reviewed by experienced radiologists, and confidentiality of any information regarding patients was maintained.
Mathematical modelling of noise in Ultrasound Image Noise in ultrasound images can be dissociated into two classes based on their effect on image pixels. One was observed to be additive and the other was observed to be multiplicative to the image pixels. Complete Mathematical modeling of additive and multiplicative noise can be represented as equation (1): Where ( , ) is the original image, corrupted with the multiplicative noise ( , ) and the additive noise ( , ).
Additive Noise removal A smoothing Gaussian filter was used for removal of the zero mean additive noise [32]. Only the application of a smoothing filter will also degrade the image quality. So here Central gradient filter was also applied to preserve the image quality. Application of Central gradient and Gaussian filters are given in equation (2) and (3) Where and are the output of Gaussian ( ) and gradients filters. The Gaussian smoothing filter ( ) of size 5×5 [33], in combination with central difference gradient kernel [34] was applied over the image as given in equation (4): Direct summation of the gradient image and smoothed image will also result in blurring. A So, sharpening factor α was introduced to avoid any blurring. the resultant image in equation (4) is free from the additive noise and only has speckle noise.

Speckle noise removal
The speckle noise abatement method proposed in this study involves multi-step processing. At first, the multiplicative noise was transformed into additive noise, by applying logarithmic transform over the enhanced image. Assuming no effect of the above process over multiplicative noise, then it can be formulated as shown in equations (5) Where ′ , are log-transformed enhanced image , desired image and the speckle noise . Then the transformed image ′ undergoes wavelet decomposition, and decomposed into four sub-bands namely , , , and . is the approximate sub-band that contains most of the important information so it remains unaltered during the rest of processing. , , and are called detail sub-band which are observed to be most affected by the speckle noise, and the rest of the processing was implemented on them. Fig-1 represents the sub-band coefficient histogram for the original image ( / / )and the noisy image ( / / ). It can be observed that the coefficient distribution of noisy mage is more stretched around 0, with lowered peak value. To recover the image the need is to inverse the transformation brought by speckle noise.

Noise Variance estimation
It is based on the fact that the speckle noise follows Rayleigh distribution [26]. The proposed estimator for speckle noise variance is a function of the scale parameter . To determine the scale parameter, the Maximum Likelihood Estimator (MLE) [35] used, as per equation (7): Where denotes coefficients of D sub-band and denotes the total number of D sub-band coefficients. The noise variance 2 can be estimated as per equation (8):

Thresholding
After estimating the noise variance next step is to define a threshold value and a thresholding function to de-speckle the image. The process for calculation is given as equation (9): Where and can be calculated from equation (10) and (11) respectively: Where, stands for the number of wavelet coefficients in the respective sub-band and is the wavelet coefficient of sub-band (D/H/V) and provides information about the decomposition level. Hereafter a three-layered distribution of sub-band coefficients and different thresholding functions was proposed as mentioned in equation (12): Where ′ is the new value assigned to coefficient of t sub-band ( / / ). The complete proposed method for noise removal is explained stepwise 1. Averaging noisy images using a Gaussian filter. 2. Generating gradient of the noisy image and summing with the averaged image as given in (11). 3. Taking log transformation of the resulting image from step 2. 4. Wavelet decomposition of logarithm transformed image, up to 2 levels. 5. Estimating noise variance from scale parameter using (15). 6. Calculate threshold value from (16). 7. Applying three-layered thresholding over each sub-band (except approximation subband). 8. Inverse wavelet decomposition of the image after thresholding. 9. At the end applying inverse logarithm transform to get a recovered image.

Results and Discussion
For performance assessment of the proposed method, the image in data sets was considered as the target image and then the noise was added to them in a controlled fashion. For evaluation of the proposed estimator of noise variance, the noise was added regarding its variance. Added noise had variance ranging from 0.01 to 0.3 with an interval of 0.01. The Sharpening factor (α) in equation (11) was set to 0.05 and the wavelet decomposition involved Daubechies waveletbased discrete wavelet transform [36]. For characterization of the performance of the proposed estimator error in estimation ( ) and percentage decrease in estimation. The error value was computed using equation (13) = | − | For evaluation of the complete denoising process other parameters like Peak Signal to Noise Ratio (PSNR), Structural Similarity Index (SSIM), and Mean Square Error (MSE) [16]were utilized. For being the highest pixel value in the original image, PSNR can be calculated as by equation (14).
For the original mage ( , ) having mean and variance , the recovered image ( , ) with mean and variance , and each image has a cross-covariance . of SSIM is calculated by the equation (15): The values of constants 1 & 2 was 6.5 and 58.5 respectively. MSE is the square of the difference of the original image and the recovered image as in equation (16).
= ∑( ( , ) − ( , )) 2 Outcomes of the proposed estimator was compared with the standard MAD estimator.  Table 2. The estimation error of the median estimator shows its lousy estimation and proves the fact that the estimator is not suited for Gaussian distribution.
The performance metrics of other practiced denoising methods named as HMBS [16], Exponential [17], Normal [19], and Wiener filter [37], were used for comparison. Fig 4 presents the PSNR value at different noise variance. For Lena image the PSNR of NormalShrink was better initially but with increasing noise variance its performance degrades rapidly. Rapid degrade in its performance with increasing noise variance was also observed in MSE and SSIM (Fig 5 and 6 respectively). While the Exponential thresholding, Weiner filter and HMBS consistently maintain their performance, but their performance also degrades very much linearly with increasing noise variance. The plots of PSNR, SSIM and MSE of proposed method tends to maintain its superiority consistently over the given range of noise variance. The recovered sub-band coefficient distribution of denoised image was shown in Fig 7, for comparison of recovery using the Soft thresholding, Hard Threshold, and Exponential thresholding. Because of the better flexibility in thresholding of the proposed multi-scale shrinkage, the recovery of coefficient histogram can be seen to be a lot improved from the other thresholding. The exponential thresholding function assigns a value to the noisy coefficient which is near zero but not zero as compares to soft thresholding where the noisy coefficients are completely grounded or assigned zero which is clear from the spike in the distribution in Fig 7. The proposed exponential shrinkage Fig 8-9 provides recovered breast ultrasound image, liver ultrasound image, Lena image, and Barbara image respectively, with different methods. Wiener filter had more effective visual recovery for Lena and Barbara images, than medical ultrasound images. The wavelet thresholding methods may have high efficiency for medical images in Fig 8 but have poor recovery for standard images in Fig 9. However, the recovery by proposed methodology has the best visual recovery among the denoising methods for Both medical and standard images.
The numerical values of PSNR, SSIM, and MSE for breast ultrasound images are given in Table 3, 4, and 5 respectively. The highest value of PSNR in Table 3 shows that the log exponential shrinkage is better at removing both speckle as well as additive noise. With the increasing noise variance, the SSIM of every method decreases rapidly except the proposed method in Table 4. It means that the processing steps in this method do not degrade the image quality as happens with the rest of the methods. The MSE of the proposed technique is lowest among them at any noise variance in Table 5. That is why the image recovered by the proposed technique is closer to the original image compared to the images recovered by HMBS, Exponential, Normal, and Wiener (Fig 8(a)-(h)). From table 3, 4, and 5 it can also be observed that the performance of the proposed method gets better from others with increasing noise variance.

Conclusion
This study implemented a wavelet thresholding method for denoising of the ultrasound image and better diagnosis. A combination of the Gaussian and gradient filters was used to remove the additive noise and preserve the details at the same time. For estimation of speckle noise, the scale parameter of the wavelet coefficient distribution was used. The threshold value was calculated using this estimated noise. For shrinkage of wavelet coefficient, a multiscale exponential function was defined. This multi scale function was more flexible in thresholding operation than any previous thresholding techniques. The complete methodology was applied to a range of medical as well as standard images. The estimation of noise variance by the designed estimator had better accuracy than the standard MAD estimator. The estimated noise variance had a very small deviation from the original variance, even at a higher noise variance. The proposed thresholding operation shows better control in coefficient thresholding than the pre-existing methods. An exponential function was used as the thresholding function. Performance evaluation of the complete methodology in terms of PSNR, SSIM, and MSE shows significant improvement from other techniques. The wavelet coefficient distribution was recovered to a large extent. The method proposed in this study is well suitable for denoising of breast ultrasound and other gray scale images. With such results, the proposed method might find application in medical image diagnosis and further empower automated medical image diagnosis.