Total variation denoising

Application of 1D total-variation denoising to a signal obtained from a single-molecule experiment.[3] Gray is the original signal, black is the denoised signal.

For a digital signal ${\displaystyle y_{n}}$, we can, for example, define the total variation as

${\displaystyle V(y)=\sum _{n}|y_{n+1}-y_{n}|.}$

Given an input signal ${\displaystyle x_{n}}$, the goal of total variation denoising is to find an approximation, call it ${\displaystyle y_{n}}$, that has smaller total variation than ${\displaystyle x_{n}}$ but is "close" to ${\displaystyle x_{n}}$. One measure of closeness is the sum of square errors:

${\displaystyle \operatorname {E} (x,y)={\frac {1}{n}}\sum _{n}(x_{n}-y_{n})^{2}.}$

So the total-variation denoising problem amounts to minimizing the following discrete functional over the signal ${\displaystyle y_{n}}$:

${\displaystyle \operatorname {E} (x,y)+\lambda V(y).}$

By differentiating this functional with respect to ${\displaystyle y_{n}}$, we can derive a corresponding Euler–Lagrange equation, that can be numerically integrated with the original signal ${\displaystyle x_{n}}$ as initial condition. This was the original approach.[1] Alternatively, since this is a convex functional, techniques from convex optimization can be used to minimize it and find the solution ${\displaystyle y_{n}}$.[3]

The regularization parameter ${\displaystyle \lambda }$ plays a critical role in the denoising process. When ${\displaystyle \lambda =0}$, there is no smoothing and the result is the same as minimizing the sum of squares. As ${\displaystyle \lambda \to \infty }$, however, the total variation term plays an increasingly strong role, which forces the result to have smaller total variation, at the expense of being less like the input (noisy) signal. Thus, the choice of regularization parameter is critical to achieving just the right amount of noise removal.

We now consider 2D signals y, such as images. The total-variation norm proposed by the 1992 article is

${\displaystyle V(y)=\sum _{i,j}{\sqrt {|y_{i+1,j}-y_{i,j}|^{2}+|y_{i,j+1}-y_{i,j}|^{2}}}}$

and is isotropic and not differentiable. A variation that is sometimes used, since it may sometimes be easier to minimize, is an anisotropic version

${\displaystyle V_{\operatorname {aniso} }(y)=\sum _{i,j}{\sqrt {|y_{i+1,j}-y_{i,j}|^{2}}}+{\sqrt {|y_{i,j+1}-y_{i,j}|^{2}}}=\sum _{i,j}|y_{i+1,j}-y_{i,j}|+|y_{i,j+1}-y_{i,j}|.}$

The standard total-variation denoising problem is still of the form

${\displaystyle \min _{y}[\operatorname {E} (x,y)+\lambda V(y)],}$

where E is the 2D L₂ norm. In contrast to the 1D case, solving this denoising is non-trivial. A recent algorithm that solves this is known as the primal dual method.[4]

Due in part to much research in compressed sensing in the mid-2000s, there are many algorithms, such as the split-Bregman method, that solve variants of this problem.

Suppose that we are given a noisy image ${\displaystyle f}$ and wish to compute a denoised image ${\displaystyle u}$ over a 2D space. ROF showed that the minimization problem we are looking to solve is:[5]

${\displaystyle \min _{u\in \operatorname {BV} (\Omega )}\;\|u\|_{\operatorname {TV} (\Omega )}+{\lambda \over 2}\int _{\Omega }(f-u)^{2}\,dx}$

where ${\textstyle \operatorname {BV} (\Omega )}$ is the set of functions with bounded variation over the domain ${\displaystyle \Omega }$, ${\textstyle \operatorname {TV} (\Omega )}$ is the total variation over the domain, and ${\textstyle \lambda }$ is a penalty term. When ${\textstyle u}$ is smooth, the total variation is equivalent to the integral of the gradient magnitude:

${\displaystyle \|u\|_{\operatorname {TV} (\Omega )}=\int _{\Omega }\|\nabla u\|\,dx}$

where ${\textstyle \|\cdot \|}$ is the Euclidean norm. Then the objective function of the minimization problem becomes:

$${\displaystyle \min _{u\in \operatorname {BV} (\Omega )}\;\int _{\Omega }\left[\|\nabla u\|+{\lambda \over 2}(f-u)^{2}\right]\,dx}$$

From this functional, the Euler-Lagrange equation for minimization – assuming no time-dependence – gives us the nonlinear elliptic partial differential equation:

For some numerical algorithms, it is preferable to instead solve the time-dependent version of the ROF equation:

$${\displaystyle {\partial u \over {\partial t}}=\nabla \cdot \left({\nabla u \over {\|\nabla u\|}}\right)+\lambda (f-u)}$$

The Rudin–Osher–Fatemi model was a pivotal component in producing the first image of a black hole.[6]

• Anisotropic diffusion
• Bounded variation
• Basis pursuit denoising
• Chambolle-Pock algorithm
• Digital image processing
• Lasso (statistics)
• Noise reduction
• Non-local means
• Signal processing
• Total variation

• TVDIP: Full-featured Matlab 1D total variation denoising implementation.
• Efficient Primal-Dual Total Variation
• TV-L1 image denoising algorithm in Matlab