A new reiterative algorithm for the Rudin-Osher-Fatemi denoising model on the graph

We consider an analogue of the well-known in image processing Rudin-Osher-Fatemi (ROF) denoising model on a general finite directed and connected graph. Then we consider the space BV of functions of bounded variation on the graph and show that the unit ball of its dual space can be described as the image of the unit ball of the space `∞ on the graph by a divergence operator. Based on this result, we propose a new fast reiterative algorithm to find the exact minimizer for the ROF model. Finally we prove convergence of the algorithm and illustrate its performance on some test examples. We would like to note that consideration of an image on the graph is important in the sense that it provides useful insights when instead of rectangular domain we have some manifold, which is a good representation for images arising in applications.


Introduction
Given a noisy observed image f ob : Ω → R in a domain Ω in R n with where f * is the original undisturbed image and η is a noise sample, the image denoising problem is how to reconstruct the best possible approximation of the original image f * from the observed image f ob .Standard variational models for image and signal denoising are based on the minimization of energy functionals consisting of a fidelity term together with a regularization term (1)(2)(3) of the type where X 0 and X 1 are Banach function spaces such that X 1 ⊂ X 0 and α, β ≥ 1.The first term in (1) is called the fidelity term and encourages minimizers that stay as close as possible to the observed image and the second term is a regularization or penalty term whose role is to try to filter out noisy oscillatory components of the observed image.The most popular choices for the function space X 1 include the space BV (Ω) of functions with bounded vari- ation and Sobolev spaces W k,2 (Ω), for k ≥ 1.The space X is usually a Lebesgue space L p (Ω) for 1 ≤ p < ∞ (see (1) ).
In this paper we consider a typical example of such functionals, known in image processing as the total variation regularization technique due to Rudin-Osher-Fatemi (ROF) (see (2)(3)(4) ) on a 2-dimensional domain Ω = [0, 1] 2 in R 2 .It suggests to take as an approximation to the original image f * the function f opt,t ∈ BV (Ω), which is the exact minimizer for the L 2,1 -functional for the couple L 2 , BV , namely for some t > 0, i.e., f opt,t ∈ BV is such that

A brief historical note on the ROF model
Since its appearance in 1992, the ROF model has received a large amount of popularity for its effeciency in regularizing images without smoothing the boundaries, and it has since been applied to a multitude of other imaging problems (see for example the book (2) ).Some of the earlier works for minimizing the total variation based on dual formulation include (5) , where the method presented is based on removal of some of the singularity caused by the non-differentiability of the quantity grad f appearing in the regularization term.In 2004, Chambolle (6) proposed an efficient projection algorithm for minimizing the total variation of an image and proved its convergence.In the last few years image decomposition models into a piecewise-smooth and oscillating components that usually people refer to as cartoon and textures + noise respectively, have received interest in the image processing community.For example f opt,t ∈ BV satisfying (3) is such that This is the decomposition of f ob into the piecewise-smooth component f opt,t ∈ BV and the component f ob − f opt,t ∈ L 2 wich contains textures and noise.The approach based on duality for such an image decomposition was introduced in 2002 by Yves Meyer in (3) by using the total variation to model the piecewise-smooth component and an appropriate dual space named G which is the Banach space composed of the distributions f which can be written f = ∂ 1 g 1 + ∂ 2 g 2 = div(g), with g 1 and g 2 in L ∞ (Ω) with the following norm to model the oscillating component.Some of the most known works proposed in the literature for numerically solving the Meyer's model or its variants include for instance the works of (7) who proposed a model that splits the image into three components, a geometrical components modeled by the total variation, a texture component modeled by a negative Sobolev norm and a noise component modeled by a negative Besov norm and the works of (8) .In 2006 (9) and (10) adapted Chambolle's algorithm in cases of presence of an operator like a convolution kernel for example.Algorithms based on the Bregman iteration are also known to be an efficient and fast way to solve total variation problems among other L 1 -regularized optimization problems.In (11) the authors designed an algorithm by using the Split Bregman iterations and the duality used by Chambolle to find the minimizer of a functional based on Meyer G-norm.Other works based on the Meyer's G-norm include for example (12) and (13) .

Preliminary notions
In this paper an analogue of (2) will be considered on the graph, and as we mentioned above, the mathematical results obtained will lead to a new algorithm.Most of the proofs of mathematical results involved will be ommitted here, interested readers can find them in (14) .The definition of BV functions can be found in (15) (see Theorem 5.3.5 on pp.227).Here we give another definition with equivalent norm.By BV we will denote the space of functions of bounded variation defined by the seminorm where is a function of y equal to variation of f on the horizontal axis for y fixed, and is a function of x equal to variation of f on the vertical axis for x fixed.
The problem of actual calculation of the function f opt,t (see ( 2) and ( 3)) is non-trivial.Standard approach is connected with discretization of the functional (2).Throughout, we will consider our discretization grid as a 2D Cartesian coordinate in (computer) screen space, i.e., with the origin at the top left corner and +y axis pointing down, which is the way matrices are represented on the computer.We divide Ω into N × N square cells and instead of the space L 2 (Ω) consider its finite dimensional subspace S N which consists of functions that are constant on each cell. where As the BV seminorm of a function g ∈ S N is equal to then the discrete analogue of the functional (2) can be written as Now we can consider the graph G = (V, E) , where the set of vertices V corresponds to cells and the set of edges E corresponds to set of pairs of cells which have common faces, i.e. and . ., N and j = 1, 2, . . ., N − 1 are the N(N − 1) horizontal edges (that we choose to be directed left right) and . ., N − 1 and j = 1, 2, . . ., N, are the N(N − 1) vertical edges (that we choose to be directed top down).So we have that dim V = N 2 and dim E = 2N (N − 1).Let us denote by S V and S E the set of real-valued functions on V and E respectively.Next we consider the analogue of the gradient operator on the graph, i.e., grad : S V −→ S E which maps function f ∈ S V to function grad f ∈ S E defined as (see Figure 1) and In this framework, the observed image f ob ∈ S N can be considered as an element of S V , and the functional (6) can be written as where , and Of course, exact minimizer of ( 7) coincides with exact minimizer of This leads to the following analogue of the ROF model on a general finite directed and connected graph G = (V, E).

Formulation of the problem on general graph
Let G = (V, E) be a finite directed and connected graph with N vertices V = {v 1 , v 2 , . . . ,v N } and M di- rected edges E = {e 1 , e 2 , . . . ,e M }, where each edge is de- termined by a pair of two vertices e = v i , v j for some i, j ∈ {1, 2, . . . ,N}.Let S V = { f : f : V −→ R} denote the set of of real-valued functions defined on the vertices and S E = {g : g : E −→ R} denote the set of of real-valued functions defined on the edges.Then the analogue of ROF model on the graph G can be formulated as follows: Problem 1. Suppose that we know function f ob ∈ S V .For given s > 0, find exact minimizer of the functional L 2,1 s, f ob ; 2 (S V ), BV(S V ) = inf g∈BV(S V ) where ; f BV(S V ) = grad f 1 (S E ) ; and (10) and operator grad : S V −→ S E is defined by the formula In this paper we will suggest a fast reiterative algorithm that constructs exact minimizer of ( 9).We would like to note that the case of a general graph could be particularly useful when instead of usually considered rectangular domain we have some manifold, for example an image on the surface of a sphere.

BV(S V )
It was shown in (16) that the exact minimizer for the L 2,1 -functional for the couple 2 , X , where space 2 is defined by the standard Euclidean norm • 2 , X is a Banach space on R n and t is a given positive parameter, is equal to the difference between f ob and the nearest element to f ob of the ball of radius t > 0 of the space X * , where X * de- notes the dual space to X.Note that proofs in (16) are also true when X is equipped with seminorm.So to construct exact minimizer for the L 2,1 -functional (9) for the couple 2 (S V ), BV(S V ) we first need to describe the ball of radius s > 0 of the space BV * (S V ) with norm defined by where Remark 1.Note that f BV(S V ) = 0 if and only if f = C for some constant C ∈ R. Therefore +∞ , otherwise.
In order to formulate the main result of this section, let us consider divergence operator on the graph, i.e. the operator div : S E → S V defined by The equation ( 13) simply means that the divergence at v is equal to the difference between the total flows into v and the total flows out of v.We can now formulate the result Theorem 3.1.The unit ball of the space BV * (S V ) is equal to the image of the unit ball of the space ∞ (S E ) under the operator div, i.e., As a corollary of Theorem 3.1 for the discrete ROF model, we will have the following result (for a reminder of space S N , see expression (4)).

Algorithm for construction of exact minimizer
Below we present a fast reiterative algorithm for the actual construction of the element f opt,t .As was mentioned in Section 3, we know from (16) that the exact minimizer f opt,t for the L 2,1 -functional is given by f opt,t = f ob − h, where h is the nearest element to f ob , in the metric of 2 (S V ), in the ball sB BV * (S V ) , i.e. inf h∈sB BV * (S V ) Next, from Theorem 3.1 we know that sB BV * (S V ) = s div B ∞ (S E ) for s > 0. The proposed algorithm constructs h through a sequence of elements g n ∈ sB ∞ (S E ) such that div(g n ) −→ h as n → +∞ in the metric of 2 (S V ).
Define the operator T as follows: where for k = 1, 2, . . ., M, T k : sB ∞ (S E ) −→ sB ∞ (S E ) is defined as follows: ,if e = e k ; g (e) ,if e = e k . where and Algorithm 1 :ROF model on the graph In the next Section we will show that div (g n ) −→ h as n → +∞ in the metric of 2 (S V ) .
To prove the convergence of the proposed algorithm we need the following proposition

Experimental results
In this section we will illustrate on some test images the performance of ROF model using our algorithm.The test observed image used has width 512 and height 512 pixels and the corresponding original image is provided for comparison purposes.The original image used (see Figure 2) is an 8-bit gray-scale image with intensity values ranging from 0 (black) to 255 (white).We notice, as expected, that the property of the ROF model for regularizing images without smoothing the boundary is preserved.
6 Below we provide another example which is a artificial image that has some geometric features with a cusp.We also notice that the noise is removed without destroying the sharp corners of the image.

Concluding remarks
In this paper, a new fast reiterative algorithm to find the exact minimizer for the ROF model is proposed and its convergence is proved.Experimental results show that the proposed algorithm is fast and does a very good job in removing the noise.It is still a remaining part of this work to study the actual rate of convergence of this algorithm and compare it with some of the traditional methods for noise reduction.

Figure 1 :
Figure 1: The neighbors of the vertex v ij .

Figure 2 :
Figure 2: The original image f * is a 512 × 512 image with intensity values ranging from 0 to 255.

Figure 4 :
Figure 4: The reconstructed image f opt,t using s = 15.5.

Figure 5 :
Figure 5: The residual image f ob − f opt,t .

Figure 6 :
Figure 6: The original image f * is a 512 × 512 image with regions of intensity values 150 and 100.

Figure 8 :
Figure 8: The reconstructed image f opt,t using s = 36.5.

Figure 9 :
Figure 9: The residual image f ob − f opt,t .