Hybrid Image Filtering Method Based on Wavelet Transform

The image acquisition process often ends with obtaining an image tainted with parasitic information which alters its quality. In order to analyze the image effectively and thus make reliable decisions, it is essential to reduce or even eliminate the noise that corrupts it. This work proposes a method of filtering images in order to reduce the noise they contain. The proposed method is a combination of bilateral filter, bivariate shrinkage and Bayes Shrink thresholding in the wavelet domain. This method consists of applying the bilateral filter to the noisy image and then decomposing it into eight wavelet sub-bands. Subsequently, bivariate shrinkage is applied to the first-level detail sub-bands, and BayesShrink thresholding is applied to the second-level detail sub-bands. The filtered image is obtained after wavelet reconstruction. The experiments were carried out on two most used images in the literature (Cameraman and Lena) corrupted with Gaussian with noise variance σ² _N = {0.05, 10, 20, 30, 50}. The performance of the proposed method was evaluated in terms of PSNR. The results obtained compared to certain approaches in the literature show that the proposed method is efficient and considerably outperforms those of existing works with an optimal PSNR value.

Introduction

An image is a signal carrying information whose analysis allows a relevant decision. Only during its acquisition can it be subject to disturbances due to environmental conditions and the equipment used. These disturbances whose role is to corrupt the nature of the image and alter its quality are none other than noise. Thus, several studies on noise have been conducted and the work of [1]–[5] allows us to define noise as a signal appearing randomly in the image and altering the quality of the latter. A noise-corrupted image does not provide a reliable interpretation of its content, so the filtering process must be used. Filtering plays a very important role regardless of the field in which it is applied because it reduces or even eliminates noise in the image. This restores the image’s quality. Several studies have been conducted to propose algorithms for image filtering. Wavelet transform-based algorithms have shown an easy determination of singularities in the image; but above all an efficient filtering of images by preserving both their structures and textures [6]–[8]. In addition, wavelet-based algorithms have resulted in better results thanks to properties such as multi-resolution and multi-scale [9]. In practice, the development of a filtering algorithm requires prior identification of the type of noises that corrupts the image. The work encountered in the literature identifies several types of noises: Gaussian noise, speckle noise, salt and pepper noise and Poisson noise. The one that draws our attention for this study is Gaussian noise. Also called electronic noise, Gaussian noise is distributed evenly over the image. It is called electronic because it appears in amplifiers or detectors. It is caused by natural sources such as the thermal vibrations of atoms and the discrete nature of radiation from hot objects [10]. In addition, it appears in an image due to factors such as electronic circuit noise and sensor noise due to poor lighting or high temperature [11]. Gaussian noise usually disturbs the gray values of digital images. This is why its model is essentially designed and characterized by its probability distribution function. This function noted $f(\tau )$ is expressed as follows: where $\tau$ represents the gray level, $n$ represents the mean of the function and $\sigma$ represents the standard deviation of the noise.

(1) $$\begin{array}{c}f\left(\tau \right)=\frac{1}{\sqrt{2\pi {\sigma }^2}}{e}^{-{(\tau -n)}^2/2{\sigma }^2}\end{array}$$

In order to obtain images of excellent quality, the problem of suppression of Gaussian noise from images has been addressed in the literature by several authors. Effective methods have been proposed. The methods resulting from these studies are generally classified into techniques based on the domain of transformation, the spatial domain, or a combination of the two methods in order to increase their benefits [11]. One of the filtering techniques in the spatial domain is the use of the bilateral filter. It is a nonlinear and local filter that has the advantage of effectively removing Gaussian noise while smoothing the image and preserving its contours [12], [13]. Although it gives better performance to eliminate noise in high-frequency areas of the image, it fails to eliminate noise in low-frequency areas. Its performance was evaluated in [12] to filter medical images. The implementation of a filtering technique based on the bilateral filter is simple but the performance of this technique depends on a number of parameters namely the window size (w) and standard deviations (${\sigma }_d$, ${\sigma }_r$). These parameters must be estimated in order to obtain an optimal result. When well chosen, the performance of the technique is better than that of linear filters such as the Wiener filter and the median filter. When the achieved performance is satisfactory, it can be used for contrast management, detail improvement, data fusion, 3D fairing, texture separation and lighting, mapping and tone management but also video processing. In order to mitigate the problems caused by Gaussian noise when acquiring a digital image, the authors of [14] proposed, applied and evaluated another filtering method. This method consists of using a threshold and histogram equalization in the wavelet domain. The performance of such a method can be improved when applying the median filter before and after the histogram equalization process. In a comparative study of the different wavelet-based filtering methods in [15], the authors showed that the wavelet transform provides an optimal representation of a noisy image. They further demonstrate that the most common way to eliminate noise using the wavelet transform is to limit the noise-affected wavelet coefficient. Thus, the threshold of the wavelet coefficient affected by noise is better only if the threshold value is correctly selected. It thus happens that the performance of wavelet-based filtering techniques depends on the estimation of the threshold value used. Thus, the advantages of wavelet transformation and thresholding were exploited in [16] to propose a filtering method that is a combination of Gaussian and bilateral filters and threshing in wavelets. This combination was implemented in order to solve the problem of loss of some image details caused by the use of the bilateral filter. Another method based on the wavelet transform was proposed in [17]. Nitin and Gupta [18] proposed a method that integrates an anisotropic Gaussian filter based on wavelets and a convolutional neural network of deep filtering. This method allowed them to eliminate noise in medical CT images. In order to develop a method that applies to character recognition, agriculture, medical science and remote sensing, Passrija et al. [19] proposed a hybridization of the Gaussian filter and the discrete wavelet transformation. This method allowed on the one hand to minimize noise in digital images and on the other hand to preserve the contours and reflect a good similarity ratio between the filtered image and the original image. Also in a hybridization approach, Rakheja and Vig [20] and Satapathy et al. [2] used a combination of the median filter and the discrete wavelet transform to filter out images corrupted by speckle and Gaussian noise. Sabah et al. [21] proposed the same combination in order to filter MRI images of the brain. Abderrahim et al. [22] suggested that the use of an embedded system would enhance the performance of an approach that combines the median filter with the wavelet transform. The work of Abderrahim et al. [22] showed that by combining the discrete wavelet transform with the fractional calculation, an efficient filtering method with better performance is obtained. Nevertheless, achieving such performance requires a good choice of fractional filter coefficients. A comparative study in [23] of different wavelet-based thresholding methods found that, for efficient filtering of images corrupted by Gaussian noise, BayesShrink’s image quality and fluidity performance are better than VisuShrink and SureShrink. Similarly, Kaur et al. [24] showed that the NeighShrinkSure threshing method is more effective than the BayesShrink method in terms of filtering Gaussian noise. In order to enhance the performance of wavelet thresholding, Fan et al. [25] and Strela et al. [26] associated the Wiener filter with wavelet threshing to filter images corrupted by Gaussian white noise. Given the fact that this is a very topical problem as researchers are working to further improve their algorithms, we propose a method for reducing Gaussian noise in images. The proposed method is a combination of bilateral filter, bivariate shrinkage and BayesShrink thresholding in the wavelet domain. This method consists of applying the bilateral filter to the noisy image and then breaking it into eight wavelet sub-bands {LL1, LH1, HL1, HH1, LL2, LH2, HL2, HH2}. Subsequently, bivariate shrinkage is applied to the first-level detail sub-bands {LH1, HL1, HH1}, and BayesShrink thresholding is applied to the second-level detail sub-bands {LH2, HL2, HH2}.

Materials and Methods

Consider $I_o$ an image without noise (original image), $N$ a Gaussian noise characterized by its variance (${{\sigma }^2}_N$), $I_b$ the resulting image of the deterioration of $I_o$ by $N$, ${\hat{I}}_o$ an approximation of the image $I_o$ (filtered image). The purpose of the proposed method is to minimize $N$ so that the system below is verified:

(2) $$\begin{array}{r}\begin{array}{c}\{\begin{array}{l}If{{\sigma }^2}_N=0then{\hat{I}}_o=I_o\\ If{{\sigma }^2}_N\ne 0then{\hat{I}}_o\to I_o\end{array}\end{array}\end{array}$$

The architecture of the proposed method which describes the filtering process is shown in Fig. 1, while the flow chart is shown in Fig. 2.

Fig. 1. Architecture of the proposed method.

Fig. 2. Flowchart of the proposed method.

After loading the original image, it is noisy. When the noise variance is zero the filtered image is identical to the original image. However, when the variance is non-zero, the bilateral filter is applied to the noisy image. The $I_{bilfilter}$ image from the bilateral filter is transformed into discrete wavelets. From this transform eight sub-band images {LL1, LH1, HL1, HH1, LL2, LH2, HL2, HH2} are obtained. Four of these sub-band images are first-level {LL1, LH1, HL1, HH1} and the four others are second-level {LL2, LH2, HL2, HH2}. After wavelet transforms, second-level details sub-bands {LH2, HL2, HH2} pass through BayesShrink thresholding while first-level {LH1, HL1, HH1} pass through bivariate shrinkage. The filtered image is obtained after applying the inverse discrete wavelet transform to the various sub-band images from the BayesShrink thresholding {$LH{2}^{{ }^{\prime }}$, $HL{2}^{{ }^{\prime }}$, $HH{2}^{{ }^{\prime }}$}, the bivariate shrinkage {$LH{1}^{{ }^{\prime }}$, $HL{1}^{{ }^{\prime }}$, $HH{1}^{{ }^{\prime }}$} and the subband image LL2. The filtering process is then evaluated using the PSNR parameter.

Bilateral Filter

The principle of the bilateral filter is to replace the intensity of each pixel with a weighted average of the intensity values of the neighboring pixels [27], [28]. In addition, it takes into account both the similarities of gray levels and the geometric proximity of neighboring pixels. Thus, pixels that look similar cause strong smoothing (homogeneous areas) while pixels that do not look similar cause weak smoothing. The behavior of the bilateral filter implemented in this work is governed by the parameters below:

$I_b$ the noisy image, ${\sigma }_S$ the standard deviation in the space domain, ${\sigma }_r$ the standard deviation in the intensity range and $w$ the window size of the filter.

The system below reflects the behavior of the implemented bilateral filter:

(3) $$\begin{array}{r}\begin{array}{c}\{\begin{array}{l}If{\sigma }_s\to \mathrm{\infty }then{\hat{I}}_{o}tendstobefuzzy\\ If{\sigma }_r\to \mathrm{\infty }then{\hat{I}}_{o}tendstoloseitscontours\end{array}\end{array}\end{array}$$

Algorithm 1 implements the bilateral filter. It is inspired by a 2D bilateral filtering function using the method described in [29] written by Lanman Douglas in [30]. This algorithm consists of 1) initializing the filter parameters $I_b$, ${\sigma }_S$, ${\sigma }_r$ and $w$. 2) Computing the Gaussian distance weights. 3) Extracting the local region. 4) Computing the Gaussian intensity weights and finally 5) Computing the bilateral filter response.

Discrete Wavelet Transform

The principle of the discrete wavelet transform is to make successive decompositions of the noisy image. Thus, the noisy image passes through a series of filter banks (a low-pass filter rated L and a high-pass filter rated H). Any decomposition is always followed by a reconstruction. Reconstruction is the inverse operation to decomposition (inverse discrete wavelet transform). In other words, a low-pass filter L and a high-pass filter H (which are conjugated filters of L and H and the same size as them) are used to perform the reconstruction. The parameters taken into account by the wavelet transform are:

$I_{bilfilter}$ the image resulting from the bilateral filter, $wm$ the mother wavelet and α the number of decompositions.

Algorithm 2 describes the discrete wavelet transform. This algorithm consists of 1) decomposing $I_{bilfilter}$ image into four wavelet sub-bands {LL1, LH1, HL1, HH1}. LH1, HL1, and HH1 are high-frequency sub-band images in horizontal, vertical, and diagonal directions, respectively, while LL1 is a low-frequency sub-band image 2) A second decomposition is carried out by applying the low-pass and high-pass filters to the image of the LL1 sub-band. This decomposition takes place using the same process as the first decomposition. After this decomposition done by dwt2 function, four sub-bands images are also obtained {LL2, LH2, HL2, HH2}. The reconstruction is carried out after thresholding and bivariate shrinkage operations.

BayesShring Thresholding

BayesShrink thresholding is a thresholding technique based on wavelets. The aim of this technique is to find a threshold value that reduces the Bayesian risk. It uses a soft threshold and depends on the wavelet sub-band, which means that thresholding is performed at each sub-band independently. It completely filters the flat regions of the image but is less sensitive to noise on the edges of the image. The Bayes threshold noted $S_{By}$ is calculated from (4):

where ${\sigma }^2$ is the noise variance and ${{\sigma }_{SB}}^2$ is the wavelet sub-band variance.

(4) $$\begin{array}{c}{S}_{By}=\frac{{\sigma }^2}{{{\sigma }_{SB}}^2}\end{array}$$

The parameters taken into account by this thresholding technique are:

$S_{By}$ the Bayes threshold value, ${\sigma }^2$ the noise variance and $SB$ the wavelet sub-band.

BayesShrink thresholding is applied to LH2, HL2 and HH2 images sub-bands from the second decomposition. Algorithm 3 describes this technique. It consists of 1) recovering the length of the sub-band image 2) estimating the variance of the sub-band ${{\sigma }_{SB}}^2.$ Then 3) computing the threshold $S_{By}$.

Bivariate Shrinkage

Bivariate shrinkage is a filtering function created in 2002 [31], [32]. Its principle is based on the use of a threshold on sub-bands of wavelets called child and parent. The parameters that govern the use of bivariate shrinkage are:

$S{B}_{child}$ the child sub-band, $S{B}_{parent}$ the parent sub-band and $T$ the threshold value.

Bivariate shrinkage is applied to LH1, HL1 and HH1 sub-band images from the first decomposition. Algorithm 4 describes this function. It consists of 1) computing the dependence $Q$ between the child sub-band and its parent, then 2) computing the shrinkage function $Q1$ and finally 3) Applying the shrinkage function to the child subbands LH1, HL1 and HH1.

The Proposed Method

Algorithm 5 is the overall algorithm of the proposed method. It is a combination of the 4 techniques previously described in 2.1, 2.2, 2.3 and 2.4 sections is summarized as follows. 1) load the original image, then 2) convert it to double and convert its pixels in the interval [0 1], then 3) generate a $noise\mathrm{\_}signal$ noisy signal from the noise variance (${{\sigma }^2}_N$). 4) Corrupt the image using $noise\mathrm{\_}signal$ noisy signal 5) Apply bilateral filter to noisy image, 6) Apply discrete wavelet transform to image $I_{bilfilter}$ resulting from the bilateral filter, 7) Compute the noise variance in level 1 sub-band LL1. 8) Compute the Bayes thresholds {$S_{By1}$, $S_{By2}$, $S_{By3}$} of the wavelet sub-bands {LH2, HL2, HH2} and limit each of these sub-bands. 9) Apply BayesShrink thresholding to level 2 sub-bands wavelets {LH2, HL2, HH2}. 10) Resize parent sub-bands {LH2, HL2, HH2}. 11) Apply bivariate shrinkage to child sub-bands {LH1, HL1, HH1}. 12) Apply the inverse discrete wavelet transform (reconstruction) by using $idwt2$ function. 13) Evaluate the performance of the filtering algorithm by calculating the PSNR parameter.

Evaluation Metric

In this work, the metric used to evaluate the performance of the proposed method in terms of Gaussian noise elimination is the Peak Signal/Noise Ratio (PSNR). It is calculated as follows:

where $M$ represents the maximum number of pixels in the image and $error$ represents the mean square error.

(5) $$\begin{array}{rl}PSNR& =20lo{g}_{10}\left(\frac{M^2}{error}\right)=20lo{g}_{10}\left(M\right)\\ & -20lo{g}_{10}\left(error\right)\end{array}$$

Work Environment

The working environment used to implement this filtering method is Matlab R2017a installed on a computer with an Intel(R)Core i7 2.90 GHz CPU and 8 GB RAM.

Results and Discussion

This section illustrates the different results obtained during the image deterioration process, the decomposition process but also the results obtained during the application of the proposed method on noisy images.

The noisy image $I_b$ is an image (cameraman or lena) corrupted using different Gaussian noise variance ${{\sigma }^2}_N$ = {0.05, 10, 20, 30, 50}. The choice of cameraman and lena images results from the fact that they belong to an open dataset and are widely used in the literature. In addition, the most recent works showing the best results in terms of image filtering (corrupted by Gaussian noise) used these images. Using these images as datasets allows us to make an excellent comparison of the experimental results obtained with those of the literature. Fig. 3 shows the original cameraman and lena images. Figs. 4 and 5 show the original cameraman and lena images and the noisy cameraman and lena images, respectively.

Fig. 3. Original images of (A) cameraman and (B) lena.

Fig. 4. Cameraman noisy image at (A) ${{\sigma }^2}_N$ = 0.05, (B) ${{\sigma }^2}_N$ = 10, (C) ${{\sigma }^2}_N$ = 20, (D) ${{\sigma }^2}_N$ = 30, and (E) ${{\sigma }^2}_N$ = 50.

Fig. 5. Lena noisy image at (A) ${{\sigma }^2}_N$ = 0.05, (B) ${{\sigma }^2}_N$ = 10 (C) ${{\sigma }^2}_N$ = 20, (D) ${{\sigma }^2}_N$ = 30, and (E) ${{\sigma }^2}_N$ = 50.

Figs. 6 and 7 represent the first and second-level sub-band images of the cameraman and lena images, respectively.

Fig. 6. Two-level decomposition process of the cameraman image (A) sub-band LL1, (B) sub-band LH1, (C) sub-band HL1, (D) sub-band HH1, (E) sub-band LL2, (F) sub-band LH2, (G) sub-band HL2, (E) sub-band HH2.

Fig. 7. Two-level decomposition process of the lena image (A) sub-band LL1, (B) sub-band LH1, (C) sub-band HL1, (D) sub-band HH1, (E) sub-band LL2, (F) sub-band LH2, (G) sub-band HL2, (E) sub-band HH2.

Tables I and II show, respectively, the results of experiments from the application of the proposed method on our dataset and the comparison of these results with some methods of the literature.

Table I. The Different Values of PSNR

Image	Noise variance ${{\mathit{\sigma }}^{\mathbf{\text{2}}}}_{\mathit{N}}$	PSNR (dB)
Cameraman	0.05	75.1974 ± 0.0001
	10	75.1980 ± 0.0064
	20	75.1209 ± 0.0159
	30	74.8470 ± 0.0506
	50	73.5457 ± 0.0343
Lena	0.05	75.7165 ± 0.0001
	10	75.7508 ± 0.0099
	20	75.8055 ± 0.0168
	30	75.7129 ± 0.0363
	50	74.5855 ± 0.0646

Table II. Comparison with Existing Works

Image	Noise variance ${{\sigma }^2}_N$			Existing methods (authors)	PSNR of existing methods (dB)	PSNR of proposed method (dB)
Cameraman	10			Discrete wavelet transform combined to median filter (Sabah et al. [21])	45.4153470	75.1980 ± 0.0064
				Haar wavelet combined to fractional filter (Abderrahim et al. [22])	24.61
	0.05			VisuShrink thresholding in wavelet domain combined to pre-gaussian filter (Nitin and Gupta [18])	22.38	75.1974 ± 0.0001
			Average		30.8018	75.1977 ± 0.0032
Lena	10			Discrete wavelet transform combined to median filter (Sabah et al. [21])	25.8820582	75.7508 ± 0.0099
				Haar wavelet combined to fractional filter (Abderrahim et al. [22])	28.43
	0.05			VisuShrink thresholding in wavelet domain combined to pre-gaussian filter (Nitin and Gupta [18])	25.80	75.7165 ± 0.0001
			Average		26.7040	75.73365 ± 0.005

Figs. 8 and 9 show the comparison between the existing methods and proposed method and the comparison results between cameraman and lena images, respectively.

Fig. 8. Comparison between the proposed method and existing methods.

Fig. 9. PSNR comparison between cameraman and lena images.

Figs. 10 and 11 show the visual rendering of cameraman and lena images after applying the proposed method.

Fig. 10. Visual rendering of cameraman image after applying the proposed method.

Fig. 11. Visual rendering of lena image after applying the proposed method.

When the noise has a variance of 10, the PSNR value is maximum for the cameraman image and for the lena image, the maximum PSNR value is obtained when the noise variance is 20. In addition, the lena image shows the best PSNR values. These results are presented, respectively in Tables I and Fig. 9. Table II shows a comparison of the proposed method with some advanced methods encountered in the literature such as the combination of the median filter and the discrete wavelet transform, the combination of the Haar wavelet with the fractional filter and the combination of the visushrink thresholding in wavelet domain and the pre-Gaussian filter. These comparisons are made on the basis of the PSNR parameter. The averages obtained show that the proposed method has approximately 44.40 (dB) more than the existing methods. In addition, Fig. 8 shows an illustration of this comparison. This comparison shows that for the different variance values used (0.05 and 10), the proposed method gives better results. Figs. 10 and 11 show the visual rendering obtained after the filtering process of the images by the proposed method. These figures each present a series of 10 images. The top 5 images are noisy images which the noise variance used is mentioned. The 5 images show the images obtained after filtering the images from above. Below these images are mentioned the PSNR values obtained after filtering. This visual rendering shows that the proposed method effectively reduces the Gaussian noise that corrupts the images.

Conclusion

In this article, we have proposed a new method of filtering images corrupted by Gaussian noise. Unlike some filtering methods, the proposed method can reduce a wide range of noise levels. ${{\sigma }^2}_N$ = {0.05, 10} for low noise levels and ${{\sigma }^2}_N$ = {20, 30, 50} for high noise levels. To evaluate the performance of our filtering method, two images widely used in the literature were exploited and the filtering performance was compared to other advanced methods in terms of PSNR. In all experimental cases (different noise variances), the proposed method gave the highest PSNR and outperformed the other methods proposed in the literature. They clearly show the effectiveness of wavelets when combined with other filtering methods. In addition, the visual quality of the images obtained makes it possible to show the effectiveness of the proposed method in terms of noise suppression but also in terms of preserving the texture and contours of the images. Given the performance of this method, it could be used in the filtering of digital images but also as a preliminary to the segmentation process.