write up for NIIII

Applied Numerical Mathematics 62 (2012) 51–66

Contents lists available at SciVerse ScienceDirect

Applied Numerical Mathematics

www.elsevier.com/locate/apnum

Correspondence between frame shrinkage and high-order nonlinear diffusion

Qingtang Jiang ∗

Department of Mathematics and Computer Science, University of Missouri–St. Louis, St. Louis, MO 63121, USA

a r t i c l e i n f o a b s t r a c t

Article history: Received 22 October 2010 Received in revised form 19 August 2011 Accepted 6 October 2011 Available online 13 October 2011

Keywords: Nonlinear diffusion filtering High-order nonlinear diffusion Signal denoising Undecimated frame filter banks Frame shrinkage Connection between nonlinear diffusion and frame shrinkage

Nonlinear diffusion filtering and wavelet/frame shrinkage are two popular methods for signal and image denoising. The relationship between these two methods has been studied recently. In this paper we investigate the correspondence between frame shrinkage and nonlinear diffusion. We show that the frame shrinkage of Ron–Shen’s continuous-linear-spline-based tight frame is associated with a fourth-order nonlinear diffusion equation. We derive high-order nonlinear diffusion equations associated with general tight frame shrinkages. These high- order nonlinear diffusion equations are different from the high-order diffusion equations studied in the literature. We also construct two sets of tight frame filter banks which result in the sixth- and eighth-order nonlinear diffusion equations. The correspondence between frame shrinkage and diffusion filtering is useful to design diffusion-inspired shrinkage functions with competitive performance. On the other hand, the study of such a correspondence leads to a new type of diffusion equations and helps to design frame-inspired diffusivity functions. The denoising results with diffusion-inspired shrinkages provided in this paper are promising.

© 2011 IMACS. Published by Elsevier B.V. All rights reserved.

1. Introduction

Nonlinear diffusion filtering [26] and wavelet shrinkage (see e.g. [15,16,22]) are two powerful methods for signal and image denoising. Correspondence between these two methods has been studied in [23,31]. In this paper we investigate the correspondence between frame shrinkage and nonlinear diffusion.

For a given 1-D signal f with a noise, nonlinear diffusion filtering is to obtain u = u(x, t) satisfying the nonlinear diffusion equation

ut = ∂

∂ x

( g ( u2x

) ux

) , (1.1)

with the initial condition

u(x, 0) = f (x), and certain boundary conditions, where g is the diffusivity and ux denotes the first-order partial derivative of u(x, t) with respect to x. The diffusivity g is a nonnegative decreasing function controlling the diffusion. The solution u(x, t) of the above nonlinear equation is a denoised version of f (x).

* Tel.: +1 314 516 6358, fax: +1 314 516 5400. E-mail address: [email protected].

0168-9274/$36.00 © 2011 IMACS. Published by Elsevier B.V. All rights reserved. doi:10.1016/j.apnum.2011.10.002

52 Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66

Since the nonlinear diffusion was introduced by Perona and Malik in 1990, a variety of nonlinear diffusion filters have been proposed, see e.g. [6,36,13] and the references therein. The fourth-order nonlinear diffusion was proposed in [38,39] to solve the problem that (the second-order) Perona and Malik diffusion and its variants tend to produce blocky effects in image denoising. The fourth-order nonlinear diffusion has also been studied in [21], and high-order diffusion with an edge enhancing functional was proposed in [35]. The theoretical properties of high-order diffusion have been studied in [14]. A 1-D high-order diffusion equation is an equation like

ut = (−1)n+1 ∂n

∂ xn

( g

(( ∂n u

∂ xn

)2) ∂n u

∂ xn

) , (1.2)

for an integer n � 2. The discretization of (1.1) could be given as follows. Let h denote the spatial step size and let τ be the time step size.

Denote

u0k = f (kh), k ∈ Z. We use u

j k , j � 1 to denote the (approximation) value of the solution u(x, t) at (kh, jτ ). Thus u

j is the approximation

solution at time jτ . With the facts that (u j+1 k − u

j k)/τ approximates ut at (kh, jτ ) and (u

j k+1 − u

j k)/h approximates ux at

(kh, jτ ), Eq. (1.1) can be discretized as

u j+1 k = u

j k +

τ

h2 g

(( u

j k+1 − u

j k

h

)2)( u

j k+1 − u

j k

) − τ h2

g

(( u

j k − u

j k−1

h

)2)( u

j k − u

j k−1

) , (1.3)

for j = 0, 1, . . .. Wavelets have been successfully used in signal and image processing [15,16,22,33]. In particular, the undecimated

wavelet transform (UWT) (also called the shift-invariant wavelet transform) based denoising [11] has been used widely for signal and image denoising. Let {p, q} be a wavelet filter bank. For a given signal {ck}k , the UWT-based denoising con- sists of the analysis step:

Ln = 1√ 2

∑ k∈Z

pk ck+n, Hn = 1√ 2

∑ k∈Z

qk ck+n, (1.4)

and the synthesis step:

uk = √

2

4

∑ n∈Z

pn Lk−n + √

2

4

∑ n∈Z

qn Sθ (Hk−n), (1.5)

where Sθ is the shrinkage function, depending a parameter θ (or several parameters). With a suitable shrinkage function (for example, the hard or soft shrinkage function), {uk}k is the denoised signal of the original signal {ck}k with noise.

It was shown in [23] that when p, q are the Haar filter pair, namely, p0 = p1 = 1, q0 = 1, q1 = −1, pk = 0, qk = 0, k �= 0, 1, then uk in (1.5) is u1k in (1.3) provided that shrinkage function Sθ and the diffusivity g satisfy

Sθ (x) = x (

1 − 4τ h2

g

( 2x2

h2

)) , (1.6)

where θ is the parameter with the diffusivity g. Namely, iterated Haar wavelet shrinking and the 2nd-order diffusion filtering result in the same signal. This relationship reveals the connection between nonlinear diffusion filtering and wavelet shrinkage and hence, it opens the gate of exchanging ideas between these two fields. In particular, the connection helps to choose shrinkage functions from diffusivity functions, and vice versa. Refer to [23,31,24] for the detailed discussion on the importance of the relationship. The reader is referred to [1,10] for the relationship between PDE diffusion and the bilateral filter, another popular method for image denoising.

Recently wavelet frames have been successfully used in noise removal [30], image recovery [7,8], image inpaint- ing/restoration [3–5], signal classification [9] and medical image analysis [18,25]. Compared with wavelet systems, the elements in a frame system may be linearly dependent; namely, frames can be redundant. The property of redundancy not only provides a flexibility for the construction of framelets with desirable properties, but also provides high sparsity of frame transform coefficients. Such sparsity is a key property for many applications. In addition, frames work better in a noisy environment [9]. It is very natural to ask whether there is a correspondence between frame shrinkage functions and the nonlinear diffusivity functions of some diffusion equations. In this paper we show that the undecimated frame shrinking corresponds to a high-order nonlinear diffusion such as

ut = ∂

∂ x

( g1

( u2x

) ux

) − ∂ 2 ∂ x2

( g2

( u2xx

) uxx

) , (1.7)

with f as initial condition:

Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66 53

u(x, 0) = f (x), where uxx denotes the second-order partial derivative of u(x, t) with respect to x. Observe that the high-order diffusion equation corresponding to a frame shrinkage is different from the high-order diffusion equations like (1.2) considered in [38,39,21,14].

The rest of the paper is organized as follows. In Section 2, we show how Ron–Shen’s continuous-linear-spline-based tight frame shrinkage corresponds to the diffusion equation given in (1.7). In Section 3, we consider the general case. We show how the vanishing moment of a highpass filter q(�) is related to the order of a nonlinear diffusion equation and derive high-order nonlinear diffusion equations associated with general tight frame shrinkages. In Section 4, we construct two sets of tight frame filter banks which result in the 6th-order and 8th-order diffusion equations. In Section 5, we provide some experiment results. We draw the conclusion in Section 6.

2. Fourth-order diffusion and tight frame shrinkage correspondence

In this section we show how Ron–Shen’s continuous-linear-spline-based tight frame shrinkage corresponds to a 4th-order nonlinear diffusion equation.

2.1. Ron–Shen’s tight frame shrinkage

For a sequence {pk}k∈Z of real numbers, we use p(ω) to denote its symbol (also called filter here):

p(ω) = 1 2

∑ k∈Z

pk e −ikω

.

Let {p, q(1), . . . , q(L)} and {p̃, q̃(1), . . . , q̃(L)} be a pair of FIR frame filter banks. Assume that they are biorthogonal, namely,

p(ω)p̃(ω) + L∑

�=1 q(�)(ω)q̃(�)(ω) = 1, (2.1)

p(ω)p̃(ω + π ) + L∑

�=1 q(�)(ω)q̃(�)(ω + π ) = 0. (2.2)

If a filter bank {p, q(1), . . . , q(L)} satisfies (2.1) and (2.2) with p̃ = p, q̃(�) = q(�), 1 � � � L, then it is called a tight frame filter bank. It was shown in [28] that if compactly supported scaling functions φ, φ̃ corresponding to lowpass filters p, p̃ are in L2(R) with

∫ R

φ(x) dx �= 0, ∫ R

φ̃(x) dx �= 0, and p(0) = p̃(0) = 1, p(π ) = p̃(π ) = q(�)(0) = q̃(�)(0) = 0, 1 � � � L, then biorthogonal frame filter banks generate wavelet bi-frames (also called dual wavelet frames) of L2(R).

Let {ck}k be the initial data. The undecimated frame transform (UFT) based denoising consists of the analysis step:

Ln = 1√ 2

∑ k∈Z

pk ck+n, H (�) n =

1√ 2

∑ k∈Z

q (�)

k ck+n, n ∈ Z, � = 1, . . . , L, (2.3)

and the synthesis step:

uk = √

2

4

∑ n∈Z

p̃n Lk−n + √

2

4

L∑ �=1

∑ n∈Z

q̃ (�) n S

� θ�

( H

(�)

k−n ) , (2.4)

where S� θ�

, 1 � � � L are the shrinkage functions, depending on parameters θ�. One can easily verify that when S �θ� (x) = x, 1 � � � L, uk is ck provided that {p, q(1), . . . , q(L)} and {p̃, q̃(1), . . . , q̃(L)} satisfy (2.1). Namely, in this case the synthesis step recovers the original signal.

In this section we consider a particular tight frame filter bank from [27]. The corresponding scaling function is the continuous linear spline function (hat function) supported on [−1, 1]. In this paper we call this filter bank Ron–Shen’s tight frame filter bank. The nonzero coefficients of the filters are

p0 = 1, p1 = p−1 = 1

2 , q

(1) 0 = 0,

q (1) −1 =

√ 2

2 , q

(1) 1 = −

√ 2

2 , q

(2) 0 = 1, q(2)−1 = q(2)1 = −

1

2 . (2.5)

With this tight frame filter bank, Ln , H (1) n , H

(2) n defined by (2.3) are

Ln = √

2 (cn−1 + 2cn + cn+1), H (1)n =

1 (cn−1 − cn+1), H (2)n =

√ 2 (2cn − cn−1 − cn+1). (2.6)

4 2 4

54 Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66

Let S 1 θ

and S 2σ denote the frame shrinkage operators applied to the first and second highpass outputs {H (1)n }n and {H (2)n }n respectively. Then the denoised signal uk after the synthesis step (2.4) is

uk = √

2

4

( Lk +

1

2 Lk−1 +

1

2 Lk+1

) +

√ 2

4

√ 2

2

( S 1θ

( H

(1) k+1

) − S 1θ (H (1)k−1))

+ √

2

4

( S 2σ

( H

(2) k

) − 1 2

S 2σ (

H (2) k+1

) − 1 2

S 2σ (

H (2) k−1

)) .

With Lk , H (1) k and H

(2) k given by (2.6), uk can be written as

uk = 1

16 (ck−2 + 4ck−1 + 6ck + 4ck+1 + ck+2) +

1

4 S 1θ

( ck − ck+2

2

) − 1

4 S 1θ

( ck−2 − ck

2

)

+ √

2

4 S 2σ

( √ 2

4 {2ck − ck−1 − ck+1}

) −

√ 2

8 S 2σ

( √ 2

4 {2ck−1 − ck−2 − ck}

)

− √

2

8 S 2σ

( √ 2

4 {2ck+1 − ck − ck+2}

) . (2.7)

With suitable shrinkage functions S 1 θ

and S 2σ , uk is the denoised signal after one step of frame denoising process of the original ck with noise. We can apply the above denoising process to uk to get further denoised signal. In fact we can apply the frame shrinkage process repeatedly to the denoised signal to get further denoised signal. We call this process the iterated frame denoising process. In the next subsection we show that the output after iterated denoising process with Ron–Shen’s tight frame filter bank is the same signal resulted by the nonlinear diffusion of a 4th-order diffusion equation.

2.2. Fourth-order nonlinear diffusion equation

We consider nonlinear diffusion equation (1.7) for u = (x, t) with u(x, 0) = f (x). To discretize the diffusion equation (1.7), we recall two formulas to approximate derivatives of a function. For a function L(x) on R and ε > 0, we have that (see e.g. [2])

L′(x0) = 1

2ε

( L(x0 + ε) − L(x0 − ε)

) − ε2 6

L(3)(ξ1), (2.8)

provided that L ∈ C 3[x0 − ε, x0 + ε], where ξ1 ∈ [x0 − ε, x0 + ε]; and that

L′′(x0) = 1

ε2

( L(x0 − ε) − 2L(x0) + L(x0 + ε)

) − ε2 12

L(4)(ξ2), (2.9)

provided that L ∈ C 4[x0 − ε, x0 + ε], where ξ2 ∈ [x0 − ε, x0 + ε]. Next we discretize (1.7) by using (2.8) and (2.9) to approximate the first and second partial derivatives with respect to

the variable x in (1.7). Recall that h and τ denote the spatial step size and the time step size respectively. As in Section 1, we use u j to denote the approximation solution of (1.7) at time jτ , and u

j k to denote the approximation value of the

solution at (kh, jτ ). Thus (u j+1 k − u

j k)/τ is the approximation of ut at (kh, jτ ). We use (2.8) with ε = h to approximate the

first-order partial derivative ux of u(x, t) and the first-order partial derivative of g1(u 2 x )ux in the first term on the right-hand

side of (1.7), while we use (2.9) with ε = h to approximate uxx and the second-order partial derivative with respect to x of g2(u

2 xx)uxx in the second term on the right-hand side of (1.7). Then we have

u j+1 k − u

j k

τ = 1

2h

{ g1

(( u

j k+2 − u

j k

2h

)2) u jk+2 − u jk 2h

− g1 ((

u j k − u

j k−2

2h

)2) u jk − u jk−2 2h

}

− 1 h2

{ g2

(( u

j k−2 − 2u

j k−1 + u

j k

h2

)2) u jk−2 − 2u jk−1 + u jk h2

− 2g2 ((

u j k−1 − 2u

j k + u

j k+1

h2

)2) u jk−1 − 2u jk + u jk+1 h2

+ g2 ((

u j k − 2u

j k+1 + u

j k+2

h2

)2) u jk − 2u jk+1 + u jk+2 h2

} .

Thus, we get

Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66 55

u j+1 k = u

j k +

τ

4h2 g1

(( u

j k+2 − u

j k

)2 / ( 4h2

))( u

j k+2 − u

j k

) − τ

4h2 g1

(( u

j k − u

j k−2

)2 / ( 4h2

))( u

j k − u

j k−2

) − τ

h4 g2

(( u

j k−2 − 2u

j k−1 + u

j k

)2 /h4

)( u

j k−2 − 2u

j k−1 + u

j k

) + 2τ

h4 g2

(( u

j k−1 − 2u

j k + u

j k+1

)2 /h4

)( u

j k−1 − 2u

j k + u

j k+1

) − τ

h4 g2

(( u

j k − 2u

j k+1 + u

j k+2

)2 /h4

)( u

j k − 2u

j k+1 + u

j k+2

) . (2.10)

Next, we obtain that uk in (2.7) after 1-step frame shrinkage is u 1 k after 1-step diffusing if S

1 θ

and S 2σ are related to g1(x) and g2(x) respectively as given in the next theorem.

Theorem 1. Let uk in (2.7) be the resulting signal after 1-step frame shrinking with input ck = f (kh), k ∈ Z and u1k in (2.10) be the signal after 1-step diffusing with the initial input u0k = f (kh), k ∈ Z. If

S 1θ (x) = x (

1 − 2τ h2

g1

( x2

h2

)) , S 2σ (x) = x

( 1 − 16τ

h4 g2

( 8x2

h4

)) , (2.11)

then uk = u1k for all k.

Proof. With u0k = ck , u1k in (2.10) after 1-step diffusion is

u1k = ck + τ

4h2 g1

( (ck+2 − ck)2/

( 4h2

)) (ck+2 − ck)

− τ 4h2

g1 ( (ck − ck−2)2/

( 4h2

)) (ck − ck−2)

− τ h4

g2 ( (ck−2 − 2ck−1 + ck)2/h4

) (ck−2 − 2ck−1 + ck)

+ 2τ h4

g2 ( (ck−1 − 2ck + ck+1)2/h4

) (ck−1 − 2ck + ck+1)

− τ h4

g2 ( (ck − 2ck+1 + ck+2)2/h4

) (ck − 2ck+1 + ck+2).

Write ck as

ck = 1

16 (ck−2 + 4ck−1 + 6ck + 4ck+1 + ck+2) −

1

8 (ck+2 − ck) +

1

8 (ck − ck−2)

+ 1 16

(ck−2 − 2ck−1 + ck) − 1

8 (ck−1 − 2ck + ck+1) +

1

16 (ck − 2ck+1 + ck+2).

Then we have that

u1k = 1

16 (ck−2 + 4ck−1 + 6ck + 4ck+1 + ck+2)

+ (

τ

4h2 g1

( (ck+2 − ck)2/

( 4h2

)) − 1 8

) (ck+2 − ck)

− (

τ

4h2 g1

( (ck − ck−2)2/

( 4h2

)) − 1 8

) (ck − ck−2)

+ (

1

16 − τ

h4 g2

( (ck−2 − 2ck−1 + ck)2/h4

)) (ck−2 − 2ck−1 + ck)

− 2 (

1

16 − τ

h4 g2

( (ck−1 − 2ck + ck+1)2/h4

)) (ck−1 − 2ck + ck+1)

+ (

1

16 − τ

h4 g2

( (ck − 2ck+1 + ck+2)2/h4

)) (ck − 2ck+1 + ck+2). (2.12)

Comparing (2.7) with (2.12), we obtain that uk in (2.7) after 1-step frame shrinking is u 1 k after 1-step diffusing if S

1 θ

, S 2σ , g1(x) and g2(x) satisfy (2.11). �

From Theorem 1, we immediately have the following corollary.

56 Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66

Corollary 1. With diffusivity functions g1(x), g2(x) and shrinkage functions S 1 θ

, S 2σ satisfying (2.11), iterated frame shrinking with Ron–Shen’s tight frame filer bank and nonlinear diffusing with (1.7) result in the same signal.

The correspondence (2.11) between frame shrinkage and diffusion filtering is useful to design diffusion-inspired shrinkage functions for frame signal denoising. On the other hand, this correspondence is useful to design frame-inspired diffusivity functions. In the following we give the corresponding shrinkage functions S 1

θ , S 2σ when diffusivity functions g1 , g2 are the

Perona–Malik diffusivity and Weickert diffusivity functions, and provide the associated diffusivity functions g1, g2 when S 1

θ , S 2σ are the hard shrinkage and soft shrinkage functions. The reader is referred to [23,31,24] for more diffusivity and

shrinkage functions. Assume the spatial step size h = 1. Corresponding to the Perona–Malik diffusivity [26]

g ( x2

) = c 1 + (x/λ)2 ,

where c is a constant, shrinkage functions S 1 θ

, S 2σ are

S 1θ (x) = x (

1 − 2τ c1 1 + (x/θ)2

) , S 2σ (x) = x

( 1 − 16τ c2

1 + (2 √

2x/σ )2

) ; (2.13)

while corresponding to the TV diffusivity [29], shrinkage functions S 1 θ

, S 2σ are

S 1θ (x) = x − 2τ sgn(x), S 2σ (x) = x − 4 √

2τ sgn(x). (2.14)

If g1, g2 are the Weickert diffusivity [36] given by

g ( x2

) = {

1, if x = 0, 1 − exp(−3.31488λ8/x8), if x �= 0,

then the corresponding shrinkage functions S 1 θ

, S 2σ are

S 1θ (x) = {

0, if x = 0, x(1 − 2τ + 2τ exp(−3.31488θ 8/x8)), if x �= 0, (2.15)

S 2σ (x) = {

0, if x = 0, x(1 − 16τ + 16τ exp(−3.31488σ 8/(2

√ 2x)8)), if x �= 0. (2.16)

If S 1 θ

, S 2σ are the hard shrinkage functions [15,22]:

S 1θ (x) = {

0, if |x| � θ, x, if |x| > θ, S

2 σ (x) =

{ 0, if |x| � σ , x, if |x| > σ ,

then the corresponding diffusivity functions g1 , g2 are

g1 ( x2

) = { 1

2τ , if |x| � θ, 0, if |x| > θ,

g2 ( x2

) = {

1 16τ , if |x| � 2

√ 2σ ,

0, if |x| > 2 √

2σ . (2.17)

When S 1 θ

, S 2σ are the soft shrinkage functions [16]:

S 1θ (x) = {

0, if |x| � θ, x − θ sgn(x), if |x| > θ, S

2 σ (x) =

{ 0, if |x| � σ , x − θ sgn(x), if |x| > σ ,

the corresponding diffusivity functions g1 , g2 are

g1 ( x2

) = {

1 2τ , if |x| � θ,

θ 2τ |x| , if |x| > θ,

g2 ( x2

) = ⎧⎨ ⎩

1 16τ , if |x| � 2

√ 2σ ,

σ

4 √

2τ |x| , if |x| > 2 √

2σ . (2.18)

Here we should point out that the diffusivity functions g1, g2 in either (2.17) or (2.18) are not differentiable. Thus their derivatives in the original equation (1.7) should be understood as the differences in (2.10), a discretized version of (1.7).

Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66 57

2.3. Minimizer of energy functional and Euler–Lagrange equations

In this subsection we show that the nonlinear diffusion equation (1.7) is related to the Euler–Lagrange equation of a variational functional.

Let

F ( u, u′, u′′

) = (u − f )2 + αΨ ((u′)2) + αΦ((u′′)2), where u′, u′′ denote the first- and second-order derivatives of u(x). Consider the energy functional

E ( u, u′, u′′

) = b∫

a

F ( u, u′, u′′

) dx,

for α > 0. A necessary condition for E to gain the minimum is that u satisfies the Euler–Lagrange equation (see [32] at p. 245)

∂ F

∂ u − d

dx

( ∂ F

∂ u′

) + d

2

dx2

( ∂ F

∂ u′′

) = 0. (2.19)

With

∂ F

∂ x = 2(u − f ), ∂ F

∂ u′ = 2αΨ ′((u′)2)u′, ∂ F

∂ u′′ = 2αΦ′((u′′)2)u′′,

we know (2.19) is

2(u − f ) − 2α d dx

( Ψ

′((u′)2)u′) + 2α d2 dx2

( Φ

′((u′′)2)u′′) = 0. Denote g1(x) = Ψ ′(x), g2(x) = Φ′(x). Then the above equation can be written as

u − f α

= d dx

( g1

(( u′

)2) u′

) − d2 dx2

( g2

(( u′′

)2) u′′

) . (2.20)

By introducing an artificial time variable t to u(x) and letting u(x, 0) = f (x), the left-hand side of Eq. (2.20) can be under- stood as the discretization to the time variable t of ∂

∂t u(x, t) with step size α, and Eq. (2.20) is a time discretization of the nonlinear diffusion equation (1.7). The reader is referred to [14] for more detailed discussions on the relationship between Euler–Lagrange equations and high-order diffusion equations.

3. High-order diffusion and undecimated frame shrinkage correspondence

In this section we consider general frame filter banks and derive the nonlinear diffusion equations associated with them. Recall that for a pair of frame filter banks {p, q(1), . . . , q(L)} and {p̃, q̃(1), . . . , q̃(L)}, Ln and H (�)n are the outputs of initial data {ck}k after analysis algorithm (2.3), and uk is the shrunk data given by (2.4) with shrinkage functions S �θ� , 1 � � � L. Observe that if the shrinking operators S �

θ� in (2.4) are the identity operator (namely, no shrinking process is applied), then uk = ck

if and only if this pair of frame filter banks satisfy (2.1). We call {p, q(1), . . . , q(L)} and {p̃, q̃(1), . . . , q̃(L)} a pair undecimated bi-frame filter banks if they satisfy (2.1). {p, q(1), . . . , q(L)} is called an undecimated tight frame filter bank if it satisfies (2.1) with p̃ = p, q̃(1) = q(1), . . . , q̃(L) = q(L). In this section we derive high-order nonlinear diffusion equations associated with undecimated frame filter banks. In Section 3.1, we obtain a proposition which rewrites uk in a formula which is closely related to a discretized version of some high-order diffusion equations. In Section 3.2, we derive the correspondence between nonlinear diffusion equations and undecimated bi-frame filter banks.

3.1. Undecimated bi-frame shrinkage

First we have the following lemma.

Lemma 1. Let {p, q(1), . . . , q(L)} and {p̃, q̃(1), . . . , q̃(L)} be a pair of undecimated bi-frame filter banks, namely they satisfy (2.1). Then ∑ m∈Z

p̃n pn+ j = 4δ( j) − L∑

�=1

∑ n∈Z

q̃ (�) n q

(�) n+ j , j ∈ Z, (3.1)

where δ( j) denotes the Kronecker delta sequence with δ( j) = 1 if j = 0, and δ( j) = 0 if j �= 0.

58 Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66

Proof. Denote q̃(0)(ω) = p̃(ω), q(0)(ω) = p(ω). From (2.1), we have (

1

2

)2 L∑ �=0

∑ m∈Z

q (�) m e

imω ∑ n∈Z

q̃ (�) n e

−inω = 1.

Using the substitution m = n + j, we have

1

4

L∑ �=0

∑ n∈Z

∑ j∈Z

q̃ (�) n q

(�) n+ j e

i jω = 1,

which is equivalent to

L∑ �=0

∑ n∈Z

q̃ (�) n q

(�) n+ j = 4δ( j), j ∈ Z.

Thus (3.1) holds. � Proposition 1. Suppose {p, q(1), . . . , q(L)} and {p̃, q̃(1), . . . , q̃(L)} are a pair of undecimated bi-frame filter banks, namely they sat- isfy (2.1). Let uk be the resulting signal given by (2.4) after 1-step frame shrinking of ck with these filter banks. Then

uk = ck + √

2

4

L∑ �=1

∑ m∈Z

q̃ (�) m

( S �θ� (x) − x

)∣∣∣∣ x=H (�)k−m

, k ∈ Z, (3.2)

where H (�) m is defined by (2.3).

Proof. By (3.1), we know the first summation in the right-hand side of Eq. (2.4) for uk is

1

2 √

2

∑ n∈Z

p̃n Lk−n = 1

4

∑ n∈Z

p̃n ∑ m∈Z

pmcm+k−n = 1

4

∑ n∈Z

∑ j∈Z

p̃n pn+ j ck+ j

= ∑ j∈Z

( δ( j) − 1

4

L∑ �=1

∑ n∈Z

q̃ (�) n q

(�) n+ j

) ck+ j = ck −

1

4

L∑ �=1

∑ j∈Z

∑ n∈Z

q̃ (�) n q

(�) n+ j ck+ j

= ck − √

2

4

L∑ �=1

∑ n∈Z

q̃ (�) n

1√ 2

∑ j∈Z

q (�) n+ j ck+ j = ck −

√ 2

4

L∑ �=1

∑ n∈Z

q̃ (�) n H

(�)

k−n.

Thus,

uk = ck − √

2

4

L∑ �=1

∑ n∈Z

q̃ (�) n H

(�)

k−n + 1√ 2

L∑ �=1

∑ n∈Z

q̃ (�) n S

� θ�

( H

(�)

k−n ) = ck +

√ 2

4

L∑ �=1

∑ m∈Z

q̃ (�) n

( S �θ� (x) − x

)∣∣∣∣ x=H (�)k−n

,

as desired. � 3.2. High-order nonlinear diffusion equation

For a (highpass) filter q(ω) = 12 ∑

k∈Z qk e −ikω , we say that it has vanishing moment order J if

∑ k∈Z

k j qk = 0, ∀ j with 0 � j < J .

The vanishing moments of analysis highpass filters imply the annihilation of discrete polynomials in the analysis step or decomposition algorithm, which results in sparse representations of input data.

Denote

C J = 1

J ! ∑ k∈Z

k J qk. (3.3)

Clearly, if q(ω) does not have vanishing moment order J + 1, then C J �= 0. Next we have a result which can be found in [37] about using a highpass filter for the approximation of the derivative of a function.

Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66 59

Lemma 2. If an FIR filter q(ω) has vanishing moment order J (not J + 1), then for a function F (x) smooth enough, 1

C J

1

ε J

∑ k∈Z

qk F (x + kε) = F ( J )(x) + o(1), (3.4)

1

C J

(−1) J ε J

∑ k∈Z

qk F (x − kε) = F ( J )(x) + o(1), (3.5)

where C J is defined by (3.3).

Proof. Using L’Hospital’s Rule repeatedly, we have

lim ε→0

1

ε J

∑ k∈Z

qk F (x + kε) = lim ε→0

1 d

dε (ε J )

d

dε

( ∑ k∈Z

qk F (x + kε) )

= lim ε→0

1

J ε J −1 ∑ k∈Z

kqk F ′ (x + kε) = lim

ε→0 1

J ddε (ε J −1)

d

dε

( ∑ k∈Z

kqk F ′ (x + kε)

)

= · · · = lim ε→0

1

J ! ∑ k∈Z

k J qk F ( J )

(x + kε) = 1 J !

( ∑ k∈Z

k J qk

) F ( J )(x).

Thus we have (3.4). (3.5) follows from (3.4) with ε replaced by −ε. � Let {p, q(1), . . . , q(L)} and {p̃, q̃(1), . . . , q̃(L)} be a pair of frame filter banks satisfying (2.1). Assume that q̃(�) and q(�)

have vanishing moment orders α� (not α� + 1) and β� (not β� + 1) respectively. Consider the following nonlinear diffusion equation for u = (x, t):

ut = L∑

�=1 (−1)1+α� ∂

α�

∂ xα�

( g�

(( ∂β� u

∂ xβ�

)2) ∂β� u

∂ xβ�

) , (3.6)

with f as initial condition:

u(x, 0) = f (x). Again, denote u0k = f (kh), and let u

j k denote the approximation to the value u(kh, jτ ) of u(x, t) at (kh, jτ ), where h and

τ are the spatial step size and the time step size. For the �-th term in (3.6), we use following formulas to approximate partial derivatives ∂

β� u ∂ xβ�

and ∂ α�

∂ xα� G(x, t), where G(x, t) := g�(( ∂β� u

∂ xβ� )2)

∂β� u ∂ xβ�

:

∂β� u

∂ xβ� (kh, jτ ) ≈ 1

Cβ�

1

hβ�

∑ n∈Z

q (�) n u(kh + nh, jτ ) ≈

1

Cβ�

1

hβ�

∑ n∈Z

q (�) n u

j n+k, (3.7)

∂α�

∂ xα� G(kh, jτ ) ≈ (−1)

α�

C̃α�

1

hα�

∑ m∈Z

q̃ (�) m G(kh − mh, jτ ), (3.8)

where Cβ� and C̃α� are the constants defined by (3.3) with q (�) and q̃(�) respectively. Observe that (3.7) and (3.8) follow

from (3.4) and (3.5) respectively with ε = h. With (3.7) and (3.8), (3.6) can be discretized as

u j+1 k = u

j k + τ

L∑ �=1

(−1)1+α� (−1) α�

C̃α�

1

hα�

∑ m∈Z

q̃ (�) m g�

(( 1

Cβ�

1

hβ�

∑ n∈Z

q (�) n u

j n+(k−m)

)2)( 1

Cβ�

1

hβ�

∑ n∈Z

q (�) n u

j n+(k−m)

) .

In particular, with ck = u0k , the above equation for j = 0 is

u1k = ck − τ L∑

�=1

1

C̃α� h α�

∑ m∈Z

q̃ (�) m g�

(( √ 2

Cβ� h β�

H (�)

k−m )2)( √

2

Cβ� h β�

H (�)

k−m )

, (3.9)

where H (�)m is defined by (2.3). Comparing (3.2) with (3.9), we have the following result.

60 Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66

Theorem 2. Let uk be the resulting signal in (2.4) after 1-step of frame shrinking of ck = f (kh), k ∈ Z with undecimated bi-frame filter banks {p, q(1), . . . , q(L)} and {p̃, q̃(1), . . . , q̃(L)} and shrinkage functions S�

θ� . Let u1k in (3.9) be the signal after 1-step diffusing defined

above for diffusion equation (3.6) with u0k = f (kh), k ∈ Z as the initial input. If

S �θ� (x) = x (

1 − 4τ C̃α� Cβ� h

α�+β� g�

( 2x2

(Cβ� ) 2h2β�

)) , 1 � � � L, (3.10)

then uk = u1k for all k.

For the tight frame filter bank, p̃ = p and q̃(�) = q(�). In this case the nonlinear diffusion equation corresponding to the tight frame shrinking is

ut = L∑

�=1 (−1)1+β� ∂

β�

∂ xβ�

( g�

(( ∂β� u

∂ xβ�

)2) ∂β� u

∂ xβ�

) , (3.11)

with f as initial condition:

u(x, 0) = f (x). For undecimated tight frame filter banks, the formulas used to discretize partial derivatives ∂

β� u ∂ xβ�

and ∂ β�

∂ xβ� G(x, t), for

G(x, t) := g�(( ∂β� u ∂ xβ�

)2) ∂β� u ∂ xβ�

are

∂β� u

∂ xβ� (kh, jτ ) ≈ 1

Cβ� h β�

∑ n∈Z

q (�) n u

j n+k,

∂β�

∂ xβ� G(kh, jτ ) ≈ (−1)

β�

Cβ� h β�

∑ m∈Z

q (�) m G(kh − mh, jτ ).

Then u1k after 1-step diffusing is

u1k = ck − τ L∑

�=1

1

Cβ� h β�

∑ m∈Z

q (�) m g�

(( √ 2

Cβ� h β�

H (�)

k−m )2)( √

2

Cβ� h β�

H (�)

k−m )

, (3.12)

where H (�)m is defined by (2.3). Comparing (3.2) with (3.12), we have the following result, which is a special case of Theo- rem 2.

Theorem 3. Let uk be the resulting signal in (2.4) after 1-step of frame shrinking of ck = f (kh), k ∈ Z with an undecimated tight frame filter bank {p, q(1), . . . , q(L)} and shrinkage functions S�

θ� . Let u1k in (3.12) be the signal after 1-step diffusing defined above for

diffusion equation (3.11) with u0k = f (kh), k ∈ Z. If

S �θ� (x) = x (

1 − 4τ (Cβ� )

2h2β� g�

( 2x2

(Cβ� ) 2h2β�

)) , 1 � � � L, (3.13)

then uk = u1k for all k.

Theorem 3 reveals the connection between nonlinear diffusion equations and general undecimated tight frame shrink- ages. As in Section 2.3, one can show that the nonlinear diffusion equation (3.11) is related to the Euler–Lagrange equation of a variational functional.

Since bi-frame (tight frame) filter banks are undecimated bi-frame (tight frame) filter banks, all results above hold true for bi-frame (tight frame) filter banks. Next, let us look at Ron–Shen’s tight frame filter bank again to illustrate the general theorem.

Example 1. Let {p, q(1), q(2)} be Ron–Shen’s tight frame filter bank defined by (2.5). Then β1 = 1, Cβ1 = − √

2; β2 = 2, Cβ2 = − 12 . Thus S 1θ1 (x) and S 2θ2 (x) in (3.13) are

S 1θ1 (x) = x − 4τ

(Cβ1 ) 2h2

g1

( 2x2

(Cβ1 ) 2h2

) x = x − 4τ

(− √

2)2h2 g1

( 2x2

(− √

2)2h2

) x = x − 2τ

h2 g1

( x2

h2

) x,

S 2θ2 (x) = x − 4τ (−1)β2 (Cβ )2h4

g2

( 2x2

(Cβ )2h4

) x = x − 4τ (−1)

2

(− 1 )2h4 g2 (

2x2

(− 1 )2h4 )

x = x − 16τ h4

g2

( 8x2

h4

) x.

2 2 2 2

Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66 61

Therefore relationship (3.13) of the diffusivity and shrinkage functions for this tight frame bank coincides with that in (2.11) with θ1 = θ , θ2 = σ .

4. More high-order nonlinear diffusion equations

In this section we construct two sets of tight frame filter banks which result in the 6th- and 8th-order nonlinear diffusion equations. In the following, denote

z = e−iω. 4.1. Sixth-order nonlinear diffusion equation

In this subsection we construct a tight frame filter bank with three highpass filters q(�), 1 � � � 3. q(�) has vanishing mo- ment order �, and it is symmetric or antisymmetric around the origin, namely, q(�)(−ω) = q(�)(ω) or q(�)(−ω) = −q(�)(ω). We consider the filters that are supported on [−2, 2]. That is the coefficient q(�)k = 0 if |k| > 2.

First let us look at q(3). If it has vanishing moment order 3, and it is antisymmetric around the origin and supported on [−2, 2], then it can be written as

q(3)(ω) = 1 2

e0

( − 1

2 z−2 + z−1 − z + 1

2 z2

) , (4.1)

for some e0 ∈ R. The formula in Lemma 2 for the 3rd derivative L(3)(x) of a function L(x) related to such a q(3)(ω) is

L(3)(x0) = 1

ε3

(−L(x0 − 2ε) + 2L(x0 − ε) − 2L(x0 + ε) + L(x0 + 2ε)) + O (ε2). (4.2) Next we consider q(1)(ω) and q(2)(ω). We choose q(1)(ω) to be the filter given by

q(1)(ω) = 1 2

c0 12

( z−2 − 8z−1 + 8z − z2), (4.3)

where c0 ∈ R. The reason for such a choice of q(1) is that the corresponding formula for the derivative L′(x) of a function L(x) is the so-called five-point formula (see e.g. [2]):

L′(x0) = 1

12ε

( L(x0 − 2ε) − 8L(x0 − ε) + 8L(x0 + ε) − L(x0 + 2ε)

) + O (ε4). (4.4) For q(2)(ω), we hope that the corresponding formula for the derivative L′′(x) is similar to the five-point formula for L′(x):

L′′(x0) = 1

12ε2 (−L(x0 − 2ε) + 16L(x0 − ε) − 30L(x0) + 16L(x0 + ε) − L(x0 + 2ε)) + O (ε4). (4.5)

q(2)(ω) corresponding to the formula (4.5) is given by

q(2)(ω) = 1 2

d0 12

(−z−2 + 16z−1 − 30 + 16z − z2), (4.6) where d0 ∈ R. Let p(ω) be the lowpass filter given

p(ω) = 1 2

( b0 z

−2 + 1 2

z−1 + 1 − 2b0 + 1

2 z + b0 z2

) , (4.7)

where b0 ∈ R. For p, q(1), q(2) and q(3) given by (4.7), (4.3), (4.6) and (4.1), we find that we are unable to choose b0 , c0 , d0 , e0 such

that p, q(1), q(2) and q(3) form a tight frame filter bank with the resulting scaling function being in L2(R). Because of this we consider q(2)(ω) given by

q(2)(ω) = 1 2

d0 a0 − 2

(−z−2 + (a0 + 2)z−1 − 2a0 − 2 + (a0 + 2)z − z2), (4.8) where d0, a0 ∈ R, a0 �= 2. The corresponding formula for the 2nd derivative L′′(x) of a function L(x) is

L′′(x0) = 1

(a0 − 2)ε2 (−L(x0 − 2ε) + (a0 + 2)L(x0 − ε) − (2a0 + 2)L(x0)

+ (a0 + 2)L(x0 + ε) − L(x0 + 2ε) ) + O (ε2). (4.9)

Then we can choose a0 , b0 , c0 , d0 , e0 such that the resulting scaling function is in L 2(R). More precisely, if

62 Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66

b0 = 3 √

2 − 14 178

, c0 = √

2 + 16b0, d0 = 7 + 104b0

6 , e0 =

√−14 − 256b0 6

, a0 = 233 + 2314b0

21 ,

then φ is W 1.20414(R) and p, q(1), q(2), q(3) form a tight frame filter bank. The corresponding nonlinear diffusion equation is

ut = ∂

∂ x

( g1

(( ∂ u

∂ x

)2) ∂ u

∂ x

) − ∂

2

∂ x2

( g2

(( ∂ 2u

∂ x2

)2) ∂ 2u

∂ x2

) + ∂

3

∂ x3

( g3

(( ∂ 3u

∂ x3

)2) ∂ 3u

∂ x3

) . (4.10)

In the above paragraph and also in the next subsection, for s > 0, W s(R) denotes the Sobolev space which consists of all functions f on R satisfying

∫ R (1 + |ω|2)s| f̂ (ω)|2 dω < ∞, where f̂ (ω) denotes the Fourier transform of f (x). The Sobolev

exponent of a compactly supported scaling function φ can be characterized by the eigenvalues of the transition operator associated with the refinement mask of φ. The reader is referred to [17,34,19] for the characterization, and to [20] for the Matlab routines about calculating the Sobolev smoothness of φ.

For q(1), q(2) and q(3) given by (4.3), (4.8) and (4.1) respectively, one can calculate directly that their vanishing moment orders β j and the corresponding Cβ j defined by (3.3) are

β1 = 1, Cβ1 = c0; β2 = 2, Cβ2 = d0; β3 = 3, Cβ3 = e0. Thus, the relationship between S �

θ� and g� in Theorem 3 is given by

S 1θ1 (x) = x (

1 − 4τ c20h

2 g1

( 2x2

c20h 2

)) , S 2θ2 (x) = x

( 1 − 4τ

d20h 4

g2

( 2x2

d20h 4

)) , S 3θ3 (x) = x

( 1 − 4τ

e20h 6

g3

( 2x2

e20h 6

)) ,

(4.11)

where h and τ are the spatial step size and the time step size respectively. In conclusion, with the relationship in (4.11) for the diffusivity and shrinkage functions, the signal resulted from iterated

denoising with p, q(1), q(2) and q(3) of the tight frame filter bank given by (4.7), (4.3), (4.8) and (4.1) respectively and that resulted from diffusion governed by Eq. (4.10) with the discretization of the 1st, 2nd and 3rd partial derivatives given by (4.4), (4.9) and (4.2) respectively are the same.

4.2. Eighth-order nonlinear diffusion equation

In this subsection, we construct a tight frame filter bank with four symmetric/antisymmetric highpass filters q(�), 1 � � � 4. q(�) has vanishing moment order �, and all the filters constructed are supported on [−2, 2]. This filter bank results in a 8th-order nonlinear diffusion equation.

Let p, q(1), q(2) and q(3) be the filters given by (4.7), (4.3), (4.6) and (4.1) respectively. We construct the 4th highpass q(4)

to have vanishing moment order 4 and to be symmetric around the origin. Then q(4) is given by

q(4)(ω) = 1 2

f 0 ( z−2 − 4z−1 + 6 − 4z + z2), (4.12)

where f 0 ∈ R. The corresponding formula for the 4th-order derivative L(4) of a function L(x) is

L(4)(x0) = 1

ε4

( L(x0 − 2ε) − 4L(x0 − ε) + 6L(x0) − 4L(x0 + ε) + L(x0 + 2ε)

) + O (ε2). (4.13) If we choose

b0 = √

7057 − 95 192

, c0 = √

2 + 16b0, d0 = √

7 + 104b0 4

, e0 = √−14 − 256b0

6 , f 0 =

589 + 624b0 54

,

then φ is in the Sobolev space W 1.19195(R), and {p, q(1), q(2), q(3)} is a tight frame filter bank. The corresponding nonlinear diffusion equation is

ut = ∂

∂ x

( g1

(( ∂ u

∂ x

)2) ∂ u

∂ x

) − ∂

2

∂ x2

( g2

(( ∂ 2u

∂ x2

)2) ∂ 2u

∂ x2

)

+ ∂ 3

∂ x3

( g3

(( ∂ 3u

∂ x3

)2) ∂ 3u

∂ x3

) − ∂

4

∂ x4

( g4

(( ∂ 4u

∂ x4

)2) ∂ 4u

∂ x4

) . (4.14)

For q(1), q(2), q(3) and q(4) given by (4.3), (4.6), (4.1) and (4.12), one can calculate directly that their vanishing moment orders β j and the corresponding Cβ j defined by (3.3) are

β1 = 1, Cβ1 = c0; β2 = 2, Cβ2 = d0; β3 = 3, Cβ3 = e0; β4 = 4, Cβ4 = f 0.

Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66 63

Fig. 1. Top-left: Original signal S 1 ; Top-right: Noised signal with SNR = 6; Bottom-left: Denoised signal with Perona–Malik shrinkage; Bottom-right: De- noised signal with Weickert shrinkage.

Thus the relationship between S � θ�

and g� in Theorem 3 is given by

S 1θ1 (x) = x (

1 − 4τ c20h

2 g1

( 2x2

c20h 2

)) , S 2θ2 (x) = x

( 1 − 4τ

d20h 4

g2

( 2x2

d20h 4

)) ,

S 3θ3 (x) = x (

1 − 4τ e20h

6 g3

( 2x2

e20h 6

)) , S 4θ4 (x) = x

( 1 − 4τ

f 20 h 8

g4

( 2x2

f 20 h 8

)) .

With such a relation among the diffusivity and shrinkage functions, the signal resulted from iterated denoising with the tight frame filter bank given by (4.7), (4.3), (4.6), (4.1) and (4.12) and that resulted from diffusion governed by Eq. (4.14) with the discretization of the 1st to 4th partial derivatives given by (4.4), (4.5), (4.2), (4.13) are the same.

The reader is referred to [12] for a tight frame filter bank consisting of 4 highpass filters q(�), 1 � � � 4 with q(�) having vanishing moment order �. Its corresponding scaling function is the C 2 cubic spline supported on [−2, 2]; and its first two highpass filters q(1), q(2) are different from these given above (they cannot result in five-point formulas (4.4)(4.5) for the 1st- and 2nd-order derivatives), while up to constants, q(3) and q(4) are the filters given in (4.1) and (4.12).

Remark 1. In this section we construct two sets of tight frame filter banks which result in the 6th- and 8th-order non- linear diffusion equations. If we consider undecimated tight frame filter banks, then we will have more flexibility for the construction which will result in smoother scaling functions. The details related to such construction are omitted here.

5. Experimental results

We carried out various experiments of signal denoising based on Ron–Shen’s tight frame filter bank with different shrink- age functions. The overall performances of the diffusion-inspired shrinking (except for the one from TV diffusivity) are comparable with hard and soft threshold denoising. Actually, they perform slightly better. Here we provide experimental results with two “toy” signals, denoted as S 1 and S 2 .

For S 1 , which is shown on the top-left of Fig. 1, five noised signals are generated by adding zero-mean Gaussian noise five times to the original signal S 1 . Each noised signal has signal-to-noise ratio (SNR) = 6. SNR is defined as

64 Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66

Table 1 Signal denoising results with different shrinkage functions.

Shrinkage S 1 θ

= P−M S 1θ = Weickert S 1θ = Hard S 1θ = Soft S 1θ = P−M S 1θ = Weickert S 1θ = TV Method S 2σ = P−M S 2σ = Weickert S 2σ = Hard S 2σ = Soft S 2σ = Weickert S 2σ = P−M S 2σ = TV SNR (for S 1 ) 18.0560 18.0643 17.4083 17.8367 18.0589 18.0539 16.2321 SNR (for S 2 ) 26.1805 25.6973 24.4256 24.4167 26.3698 25.6202 16.2660

Fig. 2. Top-left: Original signal S 2 ; Top-right: Noised signal with SNR = 16; Bottom-left: Denoised signal with Perona–Malik shrinkage; Bottom-right: Denoised signal with Weickert shrinkage.

SNR = 20(log10 |s − s̄|2 − log10 |n|2), where s is the ideal signal and s̄ is the mean of s, and n is the noise. We apply 1-level Ron–Shen’s frame shrinking iteratively 50 times to each noised signal. We provide in Table 1 the SNRs of the denoised signals with different shrinkage functions. The SNR for each case in Table 1 is the average of the SNRs of the denoised signals of the five noised signals mentioned above. When we apply the Perona–Malik (denoted as P−M) diffusivity function, we choose h = 1, τ = 14 . We choose c1 = 1 when Perona–Malik diffusivity-based S 1

θ given in (2.13) is applied to the first highpass output H (1)n while we set c2 = 18

when Perona–Malik diffusivity-based S 2σ in (2.13) is applied to the second highpass output H (2) n . We set h = 1, τ = 14 when

Weickert diffusivity-based S 1 θ

defined by (2.15) is applied to the first highpass output H (1)n , while we use h = 1, τ = 116 when Weickert diffusivity-based S 2σ in (2.16) is applied to the second highpass output H

(2) n . The parameters θ and σ are

selected such that SNRs of the denoised signals are as big as possible. The TV-diffusivity-based S 1 θ

, S 2σ defined by (2.14) are independent of the parameters. Here we choose a smaller τ with τ = 132 . From Table 1, we know for S 1 , Perona–Malik diffusivity-inspired and Weickert diffusivity-inspired shrinkages perform slightly better than hard and soft shrinkages.

The second signal S 2 we consider is shown on the top-left of Fig. 2. Again five noised signals are generated by adding zero-mean Gaussian noise five times to S 2 . In this case each noised signal has SNR = 16 and we apply Ron–Shen’s frame shrinking iteratively 100 times to each noised signal (still 50 times for TV-diffusivity-inspired shrinking). The constant c1 , c2 and τ are chosen as above. The SNRs of the denoised signals with different shrinkage functions are provided also in Table 1. Again, for each case, the SNR for S 2 in Table 1 is the average of the SNRs of the denoised signals of five noised signals. This example also shows that Perona–Malik diffusivity-inspired and Weickert diffusivity-inspired shrinkages perform better than hard and soft shrinkages.

Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66 65

6. Conclusion and future work

In this paper we establish the correspondence between frame shrinkage functions and the diffusivity functions of certain high-order nonlinear diffusion equations. We start with the frame shrinkage based on a Ron–Shen’s continuous-linear- spline-based tight frame filter bank and obtain a 4th-order nonlinear diffusion equation associated with this filter bank. After that we derive high-order nonlinear diffusion equations associated with general tight frame filter banks. These high-order nonlinear diffusion equations are different from the high-order diffusion equations studied in the literature. In addition, we construct two sets of tight frame filter banks which result in the 6th- and 8th-order nonlinear diffusion equations. We also present signal denoising experiments with various shrinkage functions including diffusivity-inspired shrinkage functions.

The study of relationship between the frame shrinkage and diffusion filtering leads to a new type of diffusion equations. The derived relationship is useful to design diffusion-inspired shrinkage functions with competitive performance. On the other hand, the relationship is helpful to design frame-inspired diffusivity functions.

In this paper we consider frame shrinkage and diffusion filtering correspondence in the 1-D case. Our future work will include the study of the correspondence between 2-D frame shrinkage and diffusion filtering, the design of competitive diffusion-inspired shrinkage functions for image denoising and the construction of tight frame filter banks which result in nonlinear diffusion equations with good performances in image noise removal.

Acknowledgements

The author would like to thank Professor Charles K. Chui for helpful discussions on nonlinear diffusion filtering. Thanks to Dr. Haihui Wang for helpful discussions on the correspondence between wavelet shrinkage and diffusion filtering. The author also thanks two anonymous referees for their valuable comments and suggestions.

References

[1] D. Barash, Fundamental relationship between bilateral filtering, adaptive smoothing, and the nonlinear diffusion equation, IEEE Trans. Pattern Analysis and Machine Intelligence 24 (2002) 844–847.

[2] R.L. Burden, J.D. Faires, Numerical Analysis, 8th ed., Thomson Books/Cole, Belmont, CA, 2005. [3] J.F. Cai, R.H. Chan, L.X. Shen, Z.W. Shen, Restoration of chopped and nodded images by framelets, SIAM J. Sci. Comput. 24 (2008) 1205–1227. [4] J.F. Cai, R.H. Chan, Z.W. Shen, A framelet-based image inpainting algorithm, Appl. Comput. Harmonic Anal. 24 (2008) 131–149. [5] J.F. Cai, S. Osher, Z.W. Shen, Split Bregman methods and frame based image restoration, Multiscale Modeling and Simulation: A SIAM Interdisciplinary

Journal 8 (2009) 337–369. [6] F. Catté, P.-L. Lions, J.-M. Morel, T. Coll, Image selective smoothing and edge detection by nonlinear diffusion, SIAM J. Numer. Anal. 29 (1992) 182–193. [7] R.H. Chan, T.F. Chan, L.X. Shen, Z.W. Shen, Wavelet algorithms for high-resolution image reconstruction, SIAM J. Sci. Comput. 24 (2003) 1408–1432. [8] R.H. Chan, S.D. Riemenschneider, L.X. Shen, Z.W. Shen, Tight frame: An efficient way for high-resolution image reconstruction, Appl. Comput. Harmonic

Anal. 17 (2004) 91–115. [9] A. Chebira, Adaptive multiresolution frame classification of biomedical images, Ph.D. Dissertation, Carnegie Mellon University, 2008.

[10] C.K. Chui, J.Z. Wang, PDE models associated with the bilateral filter, Adv. Comput. Math. 31 (2009) 131–156. [11] R.R. Coifman, D.L. Donoho, Translation-invariant de-noising, in: Wavelets and Statistics, in: Springer Lecture Notes in Statistics, vol. 103, Springer-Verlag,

New York, 1994, pp. 125–150. [12] I. Daubechies, B. Han, A. Ron, Z.W. Shen, Framelets: MRA-based construction of wavelet frames, Appl. Comput. Harmonic Anal. 14 (2003) 1–46. [13] S. Didas, Denoising and enhancement of digital images – variational methods, integrodifferential equations, and wavelets, Ph.D. Dissertation, Saarland

University, 2008, http://www.itwm.de/bv/employees/didas/dissertation.pdf. [14] S. Didas, J. Weickert, B. Burgeth, Properties of higher order nonlinear diffusion filtering, Journal of Mathematical Imaging and Vision 35 (2009) 208–226. [15] D.L. Donoho, De-noising by soft-thresholding, IEEE Trans. Inform. Theory 41 (1995) 613–627. [16] D.L. Donoho, I.M. Johnstone, Ideal spatial adaptation via wavelet shrinkage, Biometrika 81 (1994) 425–455. [17] T. Eirola, Sobolev characterization of solutions of dilation equations, SIAM J. Math. Anal. 23 (1992) 1015–1030. [18] S.D. Gertz, B.G. Bodmann, D. Vela, M. Papadakis, et al., Three-dimensional isotropic wavelets for post-acquisitional extraction of latent images of

atherosclerotic plaque components from micro-computed tomography of human coronary arteries, Academic Radiology 17 (2007) 1509–1519. [19] R.Q. Jia, Q.T. Jiang, Spectral analysis of transition operators and its applications to smoothness analysis of wavelets, SIAM J. Matrix Anal. Appl. 24 (2003)

1071–1109. [20] Q.T. Jiang, P. Oswald, Triangular

√ 3-subdivision schemes: The regular case, J. Comput. Appl. Math. 156 (2003) 47–75.

[21] M. Lysaker, A. Lundervold, X.-C. Tai, Noise removal using fourth-order partial differential equation with applications to medical magnetic resonance images in space and time, IEEE Trans. Image Process. 12 (2003) 1579–1590.

[22] S. Mallat, A Wavelet Tour of Signal Processing, 2nd ed., Academic Press, San Diego, 1999. [23] P. Mrázek, J. Weickert, G. Steidl, Correspondences between wavelet shrinkage and nonlinear diffusion, in: L.D. Griffin, M. Lillholm (Eds.), Scale-Space

Methods in Computer Vision, in: Lecture Notes in Computer Science, vol. 2695, Springer, Berlin, 2003, pp. 101–116. [24] P. Mrázek, J. Weickert, G. Steidl, Diffusion-inspired shrinkage functions and stability results for wavelet denoising, International Journal of Computer

Vision 64 (2005) 171–186. [25] M. Papadakis, B.G. Bodmann, S.K. Alexander, D. Vela, et al., Texture based tissue characterization for high-resolution CT-scans of coronary arteries,

Comm. in Numer. Methods in Engineering 25 (2009) 597–613. [26] P. Perona, J. Malik, Scale space and edge detection using anisotropic diffusion, IEEE Trans. Pattern Analysis and Machine Intelligence 12 (1990) 629–639. [27] A. Ron, Z.W. Shen, Affine systems in L2(R

d ): The analysis of the analysis operators, J. Funct. Anal. 148 (1997) 408–447. [28] A. Ron, Z.W. Shen, Affine systems in L2(R

d ), II: Dual systems, J. Fourier Anal. Appl. 3 (1997) 617–637. [29] L. Rudin, S. Osher, E. Fatemi, Nonlinear total variation based noise removal algorithms, Physica D 60 (1992) 259–268. [30] L.X. Shen, M. Papadakis, I.A. Kakadiaris, I. Konstantinidis, I. Kouri, D. Hoffman, Image denoising using a tight frame, IEEE Trans. Image Proc. 15 (2006)

1254–1263. [31] G. Steidl, J. Weickert, T. Brox, P. Mrázek, M. Welk, On the equivalence of soft wavelet shrinkage, total variation diffusion, total variation regularization,

and SIDEs, SIAM J. Numer. Anal. 42 (2004) 686–713.

66 Q. Jiang / Applied Numerical Mathematics 62 (2012) 51–66

[32] G. Strang, Introduction to Applied Math, Wellesley–Cambridge Press, Wellesley, MA, 1986. [33] G. Strang, T. Nguyen, Wavelets and Filter Banks, Wellesley–Cambridge Press, Wellesley, 1996. [34] L. Villemoes, Energy moments in time and frequency for two-scale difference equation solutions and wavelets, SIAM J. Math. Anal. 23 (1992) 1519–

1543. [35] G.W. Wei, Generalized Perona–Malik equation for image restoration, IEEE Signal Process. Lett. 6 (1999) 165–167. [36] J. Weickert, Anisotropic Diffusion in Image Processing, B.G. Teubner, Stuttgart, 1998. [37] J. Weickert, G. Steidl, P. Mrázek, M. Welk, T. Brox, Diffusion filters and wavelets: What can they learn from each other?, in: N. Paragios, Y. Chen, O.

Faugeras (Eds.), Handbook of Mathematical Models in Computer Vision, Springer, New York, 2006, pp. 3–16. [38] Y.-L. You, M. Kaveh, Image enhancement using fourth order partial differential equations, in: 32nd Asilomar Conf. Signals, Systems, Computers, vol. 2,

1998, pp. 1677–1681. [39] Y.-L. You, M. Kaveh, Fourth-order partial differential equations for noise removal, IEEE Trans. Image Process. 9 (2000) 1723–1730.