4.4. Variational denoising problem

We want to apply the framework we have studied in this chapter to the image processing problem of denoising. Let \(f\) denote a noisy image. We assume that the noise adds higher frequencies to our image. So the natural intuition is to smooth the image to remove the noise. We incorporate these assumption in our variational problem, which yields the gradient regularization and the identity as the forward operator

(4.2)\[ E(u) = {\frac{1}{2}\int_{\R^n}{\left(u(x)-f(x)\right)^2~dx}+\frac{\alpha}{2}\int_{\R^n}{|\nabla u(x)|^2~dx}}\,.\]

The energy \(E\) is strictly convex, since the data term is strictly convex and the regularizer is convex. We have seen that the energy is differentiable leading to the Euler-Lagrange-equation, which is in one to one correspndence to the unique minimizer of \(E\)

\[0 = u(x)-f(x)-\alpha\Delta u(x) \quad \forall x\in \R\,.\]

To solve the Euler-Lagrange-equation we will employ the gradient flow approach from laast section. Hereby, a convolution with the kernel

\[\begin{pmatrix} 0 & 1 & 0 \\ 1 & -4 & 1 \\ 0 & 1 & 0 \\ \end{pmatrix}\]

can be used as a discretization of the Laplace-operator. It remains open to discuss which time discretization of the gradient flow is appropiate. We have seen the explicit Euler discretization

\[u_{t+1} = u_t - \tau (u_t-f-\alpha\Delta u_t) = ((1-\tau)\mathbf{1} + \tau\alpha \Delta)u_t +\tau f \,.\]

We can use the convolution kernel to define \(A \coloneqq ((1-\tau)\mathbf{1} + \tau\alpha \Delta)\) and \(A\in\R^{n\times m}\) is of the form:\

\[\begin{pmatrix} 1-\tau-4\alpha\tau & & \tau\alpha & \\ \tau\alpha & \ddots & & \ddots\\ & \ddots& & \end{pmatrix}.\]

To show convergence of this iteration towards the solution of the Euler-Lagrange equation, we will use Banachs fixedpoint theorem for contraction maps. We start with the definition of the later.

Definition 4.3 (Contraction Map)

Let \((\mathcal{X},d)\) be a metric space. Then a map \(T:\mathcal{X}\to \mathcal{X}\) is called a contractive mapping iff \(\exists\, q<1:d\left(T(x),T(y)\right)\leq q\cdot d(x,y),~\forall x,y\in\mathcal{X}.\)

Definition 4.4 (Banachs Fixedpoint Theorem)

Let \((\mathcal{X},\|\cdot\|)\) be a Banach space and \(T:\mathcal{X}\to \mathcal{X}\) a contractive map. Then, \(\exists !\, x^*\in \mathcal{X}:x^*=T(x^*)\), called fixpoint. Furthermore, for every starting point \(x_0\in \mathcal{X}\) the sequence \(x_{n+1}=T(x_n)\) converges to \(x^*\in \mathcal{X}\), i.e. \(x_n\xrightarrow{n\to \infty}x^*.\)

Remark 4.2 (Contractive maps)

1. One can show, that a matrix \(A\) induces a contractive map if \(||A||_2:=\max_{||x||=1}{ }||Ax||_2=\sqrt{\mu_{\text{max}}}<1\), where \(\mu_{\text{max}}\) is the largest eigenvalue of \(A^*A\).

2. One can show, that \(||A||_2\leq\sqrt{||A||_1\cdot||A||_{\infty}}\) and if \(A\) is symmetric, then \(||A||_2\leq\sqrt{||A||_1^2}=||A||_1\) holds.

Applying now Banachs fixedpoint theorem to our explicit Euler discretization yields a condition on the stepsize \(\tau\). We calculate the constant of the contraction map of our update.

\[||A||_\infty =|1-\tau-4\alpha\tau|+4\alpha \tau \overset{!}{<}1 ~\implies 0\leq \tau<\frac{1}{1+4\alpha} \land 0\leq\tau< \frac{1}{0.5+4\alpha}\]

As an alternative we can discretize the gradient flow with an implicit Euler scheme.

\[\begin{aligned} u_{t+1} &=u_t-\tau(u_{t+1}-f-\alpha \Delta u_{t+1})\\ \Leftrightarrow \left((1-\tau)\operatorname{Id}-\tau\alpha \Delta\right)u_{t+1} &=u_t+\tau f \end{aligned}\]

Setting now \(A = \left((1-\tau)\operatorname{Id}-\tau\alpha \Delta\right)\) we do not have restrictions on the time step parameter \(\tau >0\) since this scheme is known to be \(A\)-stable \(\forall \alpha >0\).

4.4.1. Problem:

Our variational model (4.2) manages very well to smooth out the noise of our image. But it smoothes the whole image as well, so it smoothes every contour and edge, thus changes data of the image, which we would want to keep. We see that the gradient regularizer is not that well suited for data which has natural discontinuities, like images, since it produces a smooth function. It is well suited for signal processing, where one would expect a smooth signal.

4.4.2. Solution:

We replace the gradient regularization wiht the TV-regularizer. Here we regularize the gradient in the \(l^1\)-sense, thus allowing jumps but still smoothing small variations (noise). One can formally derive the total variation functional

(4.3)\[\min_{u\in BV(\Omega)} \frac{1}{2}\int_\Omega (u-f)^2\, dx + \alpha\int_\Omega |\nabla u|\, dx\,,\]

and receive the non-smooth 1-Laplacian

\[0 = u - f - \alpha \operatorname{div}\left(\frac{\nabla u}{|\nabla u|}\right)\,.\]

Here one needs to be careful to apply a gradient flow approach, due to the non-smoothness of the operator if \(\nabla u(x)=0\). A common approach is to add a small \(\epsilon\) value to avoid division by zero. In the next chapter we are using tools from convex analysis to do non-smooth optimization. Such an approach is applicable to functionals like Total Variation regularization.