$\DeclareMathOperator{\arccosh}{arccosh} \DeclareMathOperator*{\argmin}{arg\,min} \DeclareMathOperator{\Exp}{Exp} \newcommand{\geo}{\gamma_{\overset{\frown}{#1,#2}}} \newcommand{\geoS}{\gamma} \newcommand{\geoD}{\gamma_} \newcommand{\geoL}{\gamma(#2; #1)} \newcommand{\gradM}{\nabla_{\M}} \newcommand{\gradMComp}{\nabla_{\M,#1}} \newcommand{\Grid}{\mathcal G} \DeclareMathOperator{\Log}{Log} \newcommand{\M}{\mathcal M} \newcommand{\N}{\mathcal N} \newcommand{\mat}{\mathbf{#1}} \DeclareMathOperator{\prox}{prox} \newcommand{\PT}{\mathrm{PT}_{#1\to#2}#3} \newcommand{\R}{\mathbb R} \newcommand{\SPD}{\mathcal{P}(#1)} \DeclareMathOperator{\Tr}{Tr} \newcommand{\tT}{\mathrm{T}} \newcommand{\vect}{\mathbf{#1}}$

# The proximal map of total variation

For given data $x\in\M^{\vect{n}}$ of dimension $\vect{n}=(n_1,\ldots,n_m)$, $m\in\mathbb N$, (usually $m=1,2$) this function computes the proximal map of the TV prior in a cyclic manner.

Let $\Grid = \{1,\ldots,n_1\}\times\{1,\ldots,n_2\}\times\cdots\times\{1,\ldots,n_m\}$ denote the pixel grid the data $x$ is defined on and

the forward neighbors. Then the TV prior is given by

Since each data item $x_{\vect{k}}$ appears in $m$ terms as forward neighbor and in $m$ terms in its “own” sum, computing the proximal map is challenging.

For the cyclic proxiamal point algorithm let

and analogously $N_{\mathrm{o},j}$ for the odd indices. Then the $2m$ sums of

contains each index at most once.

This function computes the proximal maps of the TV prior in a cyclic manner using $2m$ proximal maps, $m$ for the even sets, $m$ for the odd sets.

For such a set the proximal map can be computed in parallel, since each index only appears once, all proximal maps of distances can be computed independently.

This function can be extended by carrying a different weight $\lambda_j$ for a difference in dimension $j\in\{1,\ldots,m\}$ and can also be extended to diagonal differences $\vect{k}+\vect{e}_{j_1}+\vect{e}_{j_2}$ when $\lambda_{ij}\neq 0$ is introduced, i.e. the weights $\lambda\in \mathbb R_{>0}^{m\times m}$ are given as a(n upper triangular) matrix. The standard is, to treat this matrix as a constant value given by lambda. A vector for lambda is interpreted as a diagonal matrix.

Finally the distance term can be replaced by using the optional 'DifferenceProx'. This way both a HuberTV as well as a squared TV can be easily realized with this function.

For an implementation returning a parallel evaluation of the $2m$ versions see proxParallelTV.

### Matlab Documentation

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
% proxTV(M,x,lambda) compute all proxima of TV in a cyclic manner.
% Items containing a NaN will be initalized with their first prox call.
%
% INPUT
%     M         : a manifold
%     x         : data (manifold-valued)
%   lambda      : parameter of the prox
%                   if given as a vector, an entry for each dimension,
%                   or a matrix, the diagonal for the TV terms,
%                   the offdiagonals for diagonal differences.
%
% OPTIONAL
%   'FixedMask' : binary mask the size of data items in x to indicate fixed
%   data items
%   'DifferenceProx' : (@(x1,x2,lambda)
%                                   proxAbsoluteDifference(M,x1,x2,lambda))
%                      specify a prox for the even/odd TV term proxes, i.e.
%                      switch the classical TV by Huber.
%
% OUTPUT
%     y     : result of the (iterated) proximal maps)
% ---
% Manifold-valued Image Restoration Toolbox 1.2 | R. Bergmann | 2018-02-09