Abstract

Unlike inflexible structure of soft and hard threshold function, a unified linear matrix form with flexible structure for threshold function is proposed. Based on the unified linear flexible structure threshold function, both supervised and unsupervised subband adaptive denoising frameworks are established. To determine flexible coefficients, a direct mean-square error (MSE) minimization is conducted in supervised denoising while Stein’s unbiased risk estimate as a MSE estimate is minimized in unsupervised denoising. The SURE rule requires no hypotheses or a priori knowledge about clean signals. Furthermore, we discuss conditions to obtain optimal coefficients for both supervised and unsupervised subband adaptive denoising frameworks. Applying an Odd-Term Reserving Polynomial (OTRP) function as concrete threshold function, simulations for polynomial order, denoising performance, and noise effect are conducted. Proper polynomial order and noise effect are analyzed. Both proposed methods are compared with soft and hard based denoising technologies—VisuShrink, SureShrink, MiniMaxShrink, and BayesShrink—in denoising performance simulation. Results show that the proposed approaches perform better in both MSE and signal-to-noise ratio (SNR) sense.

1. Introduction

Signals are usually corrupted by noise in capturing and transmission stages due to environment disturbance and device error. Signal denoising has become an important research topic for a long time and a wide variety of denoising methods have been proposed. Due to their effectiveness and good performance, wavelet threshold methods have become a powerful tool for denoising problems since Donoho and several others’ fundamental works. The main purpose of these methods is to estimate a wide class of functions in some smoothness spaces from their corrupted versions [1]. The different energy distribution property between smoothness spaces’ functions and noise makes these methods effective: the energy of a function in smoothness is often concentrated on few coefficients while noise is spread on all coefficients.

Of the various wavelet threshold schemes, soft and hard based threshold methods are the most popular technologies and have been theoretically verified by Donoho and Johnstone [2]. They gave a near optimal threshold value in minimax sense, where the threshold is a function of noise variance and the number of samples . This universal thresholding method is known as VisuShrink [3]. In MSE sense, better threshold methods have been proposed including SureShrink [4] and BayesShrink [5]. SureShrink gets the threshold via minimizing Stein’s unbiased risk estimate [6]. Assuming a generalized Gaussian distribution for wavelet coefficients, BayesShrink obtains the threshold through a Bayesian framework. Different a priori assumptions about the statistical distribution of wavelet coefficients, such as Gaussian scale mixture [7], Mixture of Laplacians [8], and alpha-stable [9], are also proposed, while the dependency on statistical distribution restricts their flexibility.

Soft and hard based threshold schemes suffer their own flaws due to inherent defects of soft and hard threshold functions. For soft threshold schemes, systematical biased estimation could happen, while hard threshold schemes are less biased but less sensitive to small perturbations in the data. In addition, the more important drawback is that soft and hard threshold functions do not have continuous derivatives. Various improvements had been proposed by exploring new threshold functions [1, 1017], but the nonnegative garrote-like functions [1014] are still not differentiable. Zhang [1, 15], Nasri and Nezamabadi-pour [16], and Wu et al. [17], respectively, proposed a series of threshold functions with adjustable parameters. The differentiable property makes these threshold functions suitable for gradient based minimization problems.

Most of mentioned methods above depend on a single parameter threshold, which makes those methods very sensitive to threshold value and lack of freedom. More flexible and convenient strategies had been proposed by Luisier et al. [18] and Smith et al. [19]. Both employ an idea of linear combination: the threshold function is formed by a linear combination of a set of parameters. The parameters are determined via minimizing a MSE problem. Compared to gradient based minimization, this minimization is linear and is easy to solve. In this paper, we adopt the linear combination idea and propose a unified matrix form with flexible structure for threshold function. Based on the unified flexible structure threshold function, both supervised and unsupervised subband adaptive denoising frameworks are established. In supervised denoising, a direct MSE minimization is conducted while in unsupervised denoising, we apply Stein’s unbiased risk estimate as MSE estimation. The SURE rule requires no hypotheses or a priori knowledge about clean signals. Furthermore, we discuss conditions for optimal solution for both supervised and unsupervised denoising frameworks.

Our contributions can be summarized as follows: a unified linear matrix form with flexible structure for threshold function and a concrete OTRP function are proposed; both supervised and unsupervised denoising frameworks are established; conditions for guaranteeing optimal solution for minimizing problems are discussed and provided.

The paper is organized as follows. In Section 2, both supervised and unsupervised denoising frameworks are introduced. In Section 3, frameworks employing proposed OTRP function are compared with soft and hard based denoising methods and comparison results are presented. Conclusions are made in the final.

2. Proposed Approaches

2.1. Problem Settings and Denoising Scheme

In time domain, it is assumed that the clean signal is additively corrupted by noise to produce the noisy signal in the form of where , , and are discrete time series. Additive Gaussian white noise (AGWN) is only considered, with a zero mean and a variance; that is, .

We only consider orthonormal wavelet transform (OWT), which keeps energy conservation. In OWT domain, the AGWN is still Gaussian. Relationship of wavelet coefficients in vector form iswhere , , and . is discrete wavelet transform (DWT) matrix. Each vector is stacked by both detail and approximation wavelet coefficients as similarly shown in Figure 1. Each detail is a subband. Let denote the detail subband with a vector length of . represents the total number of detail subbands and also the number of decomposition level. There holds the relationship .

We measure the quality of denoising by mean-square error (MSE) defined aswhere is the expectation operator. An estimate of MSE is given asBy considering OWT, energy conservation property guarantees that the total MSE is a weighted sum of MSE in each individual subband, which is in the form ofwhere is a weight of corresponding subband. It only needs to minimize each individual subband MSE to get minimum of total MSE. In the rest of paper, we ignore subband index and give a unified notation MSE to represent MSE in any subband.

Wavelet domain denoising generally consists of a three-stage procedure. First, perform DWT on noisy signal. Then a threshold function is applied to wavelet domain coefficients. Finally, denoised signal is obtained through an inverse DWT (IDWT). We adopt this three-stage procedure to perform our proposed denoising methods. Each subband is adaptively denoised by supervised or unsupervised denoising method as shown in Figure 2. That makes our approaches work in a subband adaptive manner. There are two input signals in Figure 2 where represents the desired signal for supervised denoising and denotes noisy signal. Corresponding uppercase letters of and denote their subbands. Both methods output an estimate for clean signal .

2.2. Proposed Threshold Function

The general problem can be formulated as a construction of threshold function that minimizes the MSE. By employing linear combination idea, an estimate of a clean signal can be represented by a linear combination of atoms; atoms mean the columns of a matrix. Involving wavelet transform, an estimate of clean wavelet coefficients can be represented by a linear combination of atoms which are decided by noisy wavelet coefficients. A general threshold function can be demonstrated aswhere is an estimate for clean wavelet coefficients. is a matrix composed of atoms and each atom is decided by a function of noisy wavelet coefficients. Dimension of an atom is decided by subband wavelet coefficients number : that is, when denoising subband . Besides, atom dimension does not have direct relationship with atom number . A more deep connection is discussed in Section 2.5. Furthermore, atom number has nothing to do with detail subband number or decomposition level . is a weight vector with each element having a weight coefficient for corresponding atom, so is a weighted sum of atoms and can be denoted as .

In this paper, the adopted concrete threshold function is OTRP function with order , denoted asThis polynomial function could fit almost any curves along with varying polynomial order. The most important advantages are its differentiability and its linear flexible structure over soft and hard threshold functions. Clearly, (7) conforms with (6) where atom is . A comparison among soft and hard threshold function and a realization of a three-order are shown in Figure 3.

2.3. Supervised Denoising

In supervised denoising, a desired signal is available; thus, a direct MSE minimization problem can be carried out. As (5) indicates, total MSE is a weighted sum of MSE in each individual subband. It only needs to minimize a subband MSE to determine its flexible parameter in (6). With prior knowledge about desired signal , a direct MSE can be written asThe optimal parameter can be obtained through on condition that is positive semidefinite and invertible. These conditions are discussed in Section 2.5. We first derive the optimal solution form asAs seen from (10), the optimal parameter is affected by matrix and desired signal . The estimates of and are obtained via an estimate of expectation operator in (4): is obtained from ; thus, coefficient is eliminated by . So, the estimate of is obtained as

2.4. Unsupervised Denoising

In practice, when a desired signal is not available, a direct MSE minimization is impossible. Constructing an estimate of MSE seems a reasonable choice. A practical approach is Stein’s unbiased risk estimate [6]. This part mainly introduces a SURE based unsupervised denoising.

2.4.1. Stein’s Unbiased Risk Estimate (SURE)

This section mainly introduces Stein’s theorem stated in [6, 20] and its tailored version.

Theorem 1. is a normal random vector with mean and identity as covariance matrix. Let be an estimate of , where functionis (weakly) differentiable. IfthenStein’s unbiased risk estimate of MSE is

The proposed structure of threshold function does not strictly satisfy conditions of in Theorem 1. Because each element of in (6) is mapping while each element in is mapping. Therefore, a tailored version of SURE for is stated in the following theorem.

Theorem 2. is a normal random vector with mean and is covariance matrix. Let be an estimate of , where is (weakly) differentiable. Ifthen The unbiased estimate of MSE is

Proof. Combining (2) and (8), the MSE is represented asThe only unknown part in the above equation is while equality , proved in [18] ensures that can be transformed to a known part. Concrete deduction formula is in the following equation:So, MSE can be rewritten asEstimate the expectation of , , , and viaIt is obvious that (19) is an unbiased estimate of MSE.
Proof is completed.

2.4.2. Coefficients Determination

From (6), we first obtain element derivative of in the form ofSo, the formula of can be deduced asThus,Combining (6) and (18) and (26), MSE can be rewritten as is obtained through on condition that is positive semidefinite and invertible. These conditions the same as in supervised denoising are discussed in Section 2.5. Here, we first derive the optimal solution form asAn estimate of is obtained via Theorem 2. From (6), we first obtain pointwise derivative of that can be denoted asThen (19) can be rewritten asThe result can be found by setting the derivative over to zero. Still positive semidefinite and invertible must be satisfied for . Estimate for optimal solution is achieved by

2.4.3. Noise Variance Estimate

As it can be seen from (29) and (33), noise variance is needed in unsupervised denoising. In practice, there are various methods that can provide estimate for AGWN variance, such as the most popular median absolute deviation (MAD) [3], SVD-based [21], and block-based [22, 23] methods. MAD implements simply and achieves accurate estimate. For simplicity, we adopt MAD to estimate noise variance. The estimate is denoted aswhere is the mean value in corresponding subband.

2.5. Optimal Solution Guarantee Discussion

Optimal solution assurance conditions for both supervised and unsupervised denoising are discussed in this section. To guarantee optimal solution of by setting to zero, the matrix must be positive semidefinite. According to the definition of positive semidefinite matrix, is positive semidefinite when for any vector . So, it is easy to verify positive semidefinite of ; that is, for any vector , there always existsThus, it can be concluded that is always positive semidefinite. Positive semidefinite property guarantees that MSE is convex while invertible of makes the optimal solution of in the form of (10), (12) and (29), (33). It must be full rank to ensure invertible , so we can deduce that must be full row or column rank because of . In practice, we do not use too many atoms to form matrix and the number of noisy coefficients in a subband is much greater than atom number . In other words, with row number far greater than column number, is a very “thin” matrix. Thus, matrix tends to be full column rank and matrix tends to be full rank.

3. Experiments and Results

3.1. Experiment Settings

Standard experiment settings for latter simulations are given. The adopted signals in simulations are “blocks,” “quadchirp,” “multitone,” and “multiband” with AGWN. In order to explain the superiority, the proposed methods are also tested by real audio signal. All those clean signals are shown in Figure 4 with each signal having length of 8192 for the first four signals and a length of 65536 for audio signal. The adopted orthogonal wavelet is “sym8” and decomposition level is . So the detail subband has a vector length of , where can be 8192 or 65536.

In supervised denoising, coefficient is determined by matrix and desired subband in each individual subband. Desired subbands come from desired signal by OWT. Desired signal can be clean or noisy. In practice, it may only have access to noisy signal with known noise level; thus, we adopt it as desired signal to determine coefficient and then use to derive the denoised version of noisy signal with unknown noise level. A remarkable notice is that the desired signal and noisy signal are, respectively, the same clean signal corrupted by known and unknown noise level. In unsupervised denoising, coefficient is determined by matrix , noisy signal subband , and noise variance. There are no requirements for clean signal or desired signal. For simplicity, an estimate for noise variance is provided by MAD.

Applying the proposed OTRP function as concrete threshold function, supervised and unsupervised denoising methods are compared with several other available techniques, that is, VisuShrink, SureShrink, MiniMaxShrink, and BayesShrink. Wavelab 850 is applied to perform the other four techniques; their noise variance estimates are also provided by MAD. Denoising quality is measured by both MSE and SNR, which are defined as

3.2. Polynomial Order Simulation

Before carrying out the final denoising, proper polynomial order needs be explored. Due to the fact that polynomial order is equal to atom number of matrix , it needs to select proper polynomial order to guarantee full column rank of and thus full rank of . Although it has already been discussed in Section 2.5 that subband wavelet coefficients number is far greater than polynomial order and tends to be full rank, simulations are still needed. Whether is full rank or not is determined by polynomial order and subband wavelet coefficients. Thus, a probability statistical table is made through simulations. The achieved probabilities make up Table 1 where PO represents polynomial order and SB denotes subband. The probabilities reveal how often is full rank in different PO for different signals corrupted with different noise level.

Six orders of PO are tested in Table 1. Due to 5-level wavelet decomposition, there are five detail SBs to be processed. Each SB corresponds to a and different PO corresponds to different dimension of . So, there are 30 combinations shown in Table 1 header where each PO contains five SBs. In -axis direction, Table 1 reflects how PO influences matrix full rank probabilities. In -axis direction of Table 1, four signals are tested with each corrupted by different noise variance . This reflects how noise level influences matrix full rank probabilities. For example, the italic number shown in Table 1 means matrix full rank probability is . This matrix is determined by PO2 and “blocks” signal’s SB3 in noise level. All probabilities are obtained by counting matrix full rank times in a 100-trail.

From Table 1, PO1–PO3 always ensure matrix full rank for four signals in all tested noise level. Saying matrix full rank for a signal we mean matrix full rank for every subband of this signal. PO4 only cannot ensure “blocks” due to almost probability for SB5 (bold underlined numbers, bold numbers indicate the probabilities under ). For PO5 and PO6, there are more and more bold numbers than PO4. In horizontal direction, the trend is higher PO lower probability for satisfying matrix full rank. In vertical direction, higher noise level tends to lead lower probability for full rank, but PO is a dominant factor. From statistical analysis, we restrict our option from PO1 to PO3.

Visual comparisons of PO1 to PO3 in supervised denoising are shown in Figure 5. As it can be seen, in both MSE and SNR sense PO1 is inferior to PO2 and PO3 in all four signals. PO2 is inferior to PO3 in blocks and almost equals PO3 in the three remaining signals. In unsupervised denoising, similar conclusions can be made and are not repetitive state. In the end, we choose PO2 for final denoising simulation.

3.3. Denoising Performance Comparison

As mentioned above, comparisons are made among supervised and unsupervised denoising methods and VisuShrink, SureShrink, MiniMaxShrink, and BayesShrink. In both MSE and SNR sense, denoising performance comparisons in global and local views with different noise variance are presented in Figure 6. changes in range of 0.5–2.0 with step 0.1 and each data point plotted in Figure 6 is a 50-average result. A remarkable notice is that we adopt a clean signal corrupted by noise in level as desired signal for supervised denoising. Moreover, MSEs in Figure 6 are plotted in log scale and SNRs in Figure 6 are plotted in linear scale.

From global and local views of Figure 6, supervised denoising method achieves better MSE and SNR quality than unsupervised denoising method and BayesShrink for all four noisy signals. All these three methods have better MSE and SNR performance than VisuShrink, SureShrink, and MiniMaxShrink. In noisy “blocks” processing, BayesShrink is slightly better than unsupervised denoising method when is greater than 0.9 while poorer when is less than 0.9 in both MSE and SNR sense. However, unsupervised denoising method has better MSE and SNR performance than BayesShrink in noisy quadchirp, multitone, and multiband processing. Thus, a conclusion can be made that the proposed supervised and unsupervised denoising methods show their advantages.

Our approaches are also validated for real-world signals. We test supervised and unsupervised denoising approaches by an AGWN corrupted real audio signal. In both MSE and SNR sense, comparisons with changing in range of 0.1–0.5 are shown in Figure 7 and each data point plotted is still a 50-average result. We adopt clean audio signal corrupted by noise in level as desired signal for supervised denoising. In addition, MSEs are plotted in log scale and SNRs are plotted in linear scale. From global and local views of Figure 7, the proposed supervised and unsupervised denoising methods have smaller MSEs and higher SNRs than VisuShrink, SureShrink, MiniMaxShrink, and BayesShrink. A conclusion can still be made that the proposed supervised and unsupervised denoising methods show their advantages in real-world condition.

Remark. All denoising methods except for supervised denoising involved in Figures 6 and 7 require noise variance , and estimates of it are achieved by MAD for all those methods. So, proposed unsupervised method and the other four techniques work under the same condition, which reflects superiority of the proposed method. In practice, MAD with simple implementation and estimate accuracy do satisfy our requirements.

3.4. Noise Effect Analysis

In this section, different noise levels of desired signal for supervised denoising are analyzed. It only aims at real audio signal. In simulation, desired signal noise levels are set in range to with equal step . Both MSE and SNR results are shown in Figure 8 and each data point plotted in Figure 8 is averaged by 50 times. It must be noticed that on -axis represents noise level of noisy signal to be processed. All desired signal noise levels are compared with BayesShrink denoising mothed. From Figure 8 it can be seen that all to desired signal noise levels have better denoising performance than BayesShrink in MSE and SNR sense when on -axis. While when on -axis BayesShrink does better than and and worse than to in MSE and SNR sense. In addition, denoising performance becomes worse when desired signal noise level gets higher. In conclusion, if some noisy desired signal (even higher noise level than the signal to be processed) is available at hand, it is still possible to achieve good denoising performance.

4. Conclusions

In order to utilize the effectiveness of wavelet domain denoising, two subband adaptive denoising schemes were proposed in this paper. For any subband, a unified linear matrix form with flexible structure of threshold function was proposed and OTRP function was proposed for a concrete realization. Based on the unified linear flexible structure threshold function, both supervised and unsupervised subband adaptive denoising frameworks were established. A direct MSE minimization was conducted in supervised denoising while Stein’s unbiased risk estimate as MSE estimate was minimized in unsupervised denoising for flexible coefficients determination. Conditions to obtain optimal coefficients were further discussed.

Applying the concrete OTRP function, simulations for polynomial order, denoising performance, and noise effect were conducted. A full rank probability statistical table was generated in polynomial order simulation. The table reflects that PO1–PO3 always ensure full rank of . By choosing PO2, denoising performance was compared among proposed supervised and unsupervised denoising methods and VisuShrink, SureShrink, MiniMaxShrink, and BayesShrink. Comparison results demonstrate that the proposed methods have advantages over those soft and hard threshold function based denoising schemes in both simulation and real-world condition. Final simulation of noise effect analysis for supervised denoising demonstrates that supervised denoising method is robust for different desired signal noise level.

Competing Interests

The authors declare that there is no conflict of interests.

Acknowledgments

This work is fully supported by National Science and Technology Major Project of the Ministry of Science and Technology of China (no. Y6D8020801).