This algorithm computes the cyclic or randomized proximal point algorithm on a manifold for given input data , where denotes the data domain —usually a signal or pixel grid, but necessarily on some array form—, a set of proximal maps, and a stopping criterion.
Given a function this algorithm aims to find a minimizer
so in a sense it is similar to a gradient descent method, but it is also applicable if is not differentiable and there is no need to determine a step size; choosing the described in the following determine merely the speed of convergence.
If we know for , that there exists a decomposition of into atoms
where each of the is assumed to have a proximal map, that is either known in closed form or otherwise easily (approximately to be computed, i.e.
see [1] for an introduction in Euclidead space . The cyclic proximal point (CPP) algorithm was introduced on manifolds in [2]. It consists of applying the proximal maps iteratively
where we denote the inner (one cycle) iterations by a fractional index.
Parameters
While the M
denotes a manifold
class and x
denotes the
initial data both remaining parameters are functionals. The easier one is the stoppingCriterion = @(x,xold,lambda,iter)
with the for parameters , , and the current iteration .
The itertion is running while the functional returns false, i.e. the stoppingCriterion is not fulfilled and stops after a cycle if it returns true.
The easiest is of course just depending on the iteration: @(x,xold,lambda,iter) iter>1999
will stop after 2000 iterations.
The varibale proximalMaps
is a cell array containing the
proximal maps as anonymous functions. More precisely each entry of proximalMaps
is an anonymous function
1
fi = @(x,lambda) ...
evaluating the corresponding proximal map of . For a list, see below
Optional Parameters
lambdaIterate
 a functional
@(iter)
returning a value for lambda
 easiest way to set the lambda iterate: defines the
lambdaIterate
as Order
 if not specified, the above cyclic order is used. Other possibilities are
'permute'
which permutes the order of the proximal maps in each cycle randomly or'random'
which choses one proximal map at random per iteration to be evaluated. Debug
 If set to an anonymous functional depending on
@(x,xold,iter)
this is called in each cycle, for example to produce text output or save iteration data to a file. Record
 If set to a functional
@(x,iter)
returning a column vector, and calling thecyclicProximalPoint
with a second return value, this function records these column vectors in a matrix returned in the second return valuerecData
.
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
24
25
26
27
28
29
30
31
32
33
34
35
36
% CyclicProximalPoint(M,f,proximalMaps,stoppingCriterion)
% Compute the cyclic proximal point algorithm starting with data f,
% cyclic through the proximal points
% Input
% M : a manifold
% x : given (initial) data
% proximalMaps : cell arrays of functionals computing the proximal
% maps, i.e. can be called as prox = @(x,lambda)
% stoppingCriterion : a functional @(x,xold,lambda,iter) as a stopping
% criterion based on x, the last iterate xold,
% lambda,
% that indicates, then to stop, i.e., is 1 for stops.
%
% OPTIONAL
% 'lambdaIterate' : (@ (iter)) a functional @(iter) computing the
% iterates value of lambda
% 'lambda' : if no iterate (see above) is given, this value is used
% to set the iterate to @(iter) lambda/iter
% 'Order' : ('') the standard order is ascending in proximalMaps,
% 'random' takes a random proximal map each iteration
% (note that then one iteration is only _one_
% proximal map evaluated)
% 'permute' keeps a cyclic manner but permutes the
% order for each cycle anew.
% Debug : ([]) a function @ (x,xold,iter) producing debug
% output
% Record : (@(x,iter)) a function returning a column vector, if
% there's a second return value and this function is
% given, data is recorded in an array and returned
%
% OUTPUT
% y : result of the CPP algorithm
% recData : recorded data array of the by the recordfunction
% 
% ManifoldValued Image Restoration Toolbox 1.0
% R. Bergmann ~ 20180122  20180304
See also
used in
 Tutorial on how to compute a median on a manifold
 The CPP for additive first and second order TV
 The CPP for Huber TV
 an easy soptting criterion creator
References

Bertsekas, D P (2010). Incremental Gradient, Subgradient, and Proximal Methods for Convex Optimization: a Survey. Laboratory for Information and Decision Systems, MIT, Cambridge, MA
We survey incremental methods for minimizing a sum \(\sum_{i=1}^m f_i(x)\)consisting of a large number of convex component functions fi. Our methods consist of iterations applied to single components, and have proved very effective in practice. We introduce a unified algorithmic framework for a variety of such methods, some involving gradient and subgradient iterations, which are known, and some involving combinations of subgradient and proximal methods, which are new and offer greater flexibility in exploiting the special structure of fi. We provide an analysis of the convergence and rate of convergence properties of these methods, including the advantages offered by randomization in the selection of components. We also survey applications in inference/machine learning, signal processing, and largescale and distributed optimization.@techreport{Ber10, author = {Bertsekas, Dimitri P}, title = {Incremental Gradient, Subgradient, and Proximal Methods for Convex Optimization: a Survey}, year = {2010}, institution = {Laboratory for Information and Decision Systems, MIT}, address = {Cambridge, MA}, number = {LIDSP2848}, eprint = {1507.01030}, eprinttype = {arXiv} }

Bačák, M (2014). Computing medians and means in Hadamard Spaces. SIAM Journal on Optimization. 24 1542–66
The geometric median as well as the Fréchet mean of points in a Hadamard space are important in both theory and applications. Surprisingly, no algorithms for their computation are hitherto known. To address this issue, we use a splitting version of the proximal point algorithm for minimizing a sum of convex functions and prove that this algorithm produces a sequence converging to a minimizer of the objective function, which extends a recent result of Bertsekas [Math. Program., 129 (2011), pp. 163–195] into Hadamard spaces. The method is quite robust, and not only does it yield algorithms for the median and the mean, but also it applies to various other optimization problems. We, moreover, show that another algorithm for computing the Fréchet mean can be derived from the law of large numbers due to Sturm [Ann. Probab., 30 (2002), pp. 1195–1222]. In applications, computing medians and means is probably most needed in tree space, which is an instance of a Hadamard space, invented by Billera, Holmes, and Vogtmann [Adv. in Appl. Math., 27 (2001), pp. 733–767] as a tool for averaging phylogenetic trees. Since there now exists a polynomialtime algorithm for computing geodesics in tree space due to Owen and Provan [IEEE/ACM Trans. Comput. Biol. Bioinform., 8 (2011), pp. 2–13], we obtain efficient algorithms for computing medians and means of trees, which can be directly used in practice.@article{Bac14, author = {Ba{\v{c}}{\'a}k, Miroslav}, title = {Computing medians and means in {H}adamard Spaces}, sjournal = {SIAM J. Optim.}, journal = {SIAM Journal on Optimization}, volume = {24}, year = {2014}, number = {3}, pages = {15421566}, eprint = {1210.2145}, eprinttype = {arXiv}, doi = {10.1137/140953393} }