Convergence rate analysis of primal-dual splitting schemes

Introduction

Primal-dual algorithms are abstract splitting schemes that solve monotone inclusion and convex optimization problems. These schemes fully decompose problems built from sums, linear compositions, parallel sums, and infimal convolutions of simple functions so that each simple term is processed individually. This decomposition is achieved by cleverly combining primal and dual pair problems into a single inclusion problem, to which standard operator splitting algorithms can be applied. This process gives rise to algorithms that are inherently parallel or distributed and in which expensive matrix inversions can be avoided. The characteristics of primal-dual algorithms are especially desirable for large-scale applications in machine learning, image processing, distributed optimization, and control.

Primal-dual methods have a long history with many contributors, and an attempt to summarize and relate all of the contributions is beyond the scope of this paper. In this paper, we are mainly concerned with the line of work that began in and the many generalizations and enhancements of the basic framework that followed . Thus, we consider the following prototypical convex optimization problem as our guiding example:

where $\square$ denotes the infimal convolution operation (see Section 1.2), $n\in{\mathbf{N}}$, $n\geq 1$, ${\mathcal{H}}_{i}$ are Hilbert spaces for $i=0,\ldots,n$, the functions $f,g:{\mathcal{H}}_{0}\rightarrow(-\infty,\infty]$ and $h_{i},l_{i}:{\mathcal{H}}_{i}\rightarrow(-\infty,\infty]$ are closed, proper, and convex for $i=1,\cdots,n$, and $B_{i}:{\mathcal{H}}_{0}\rightarrow{\mathcal{H}}_{i}$ is a bounded linear map for $i=1,\ldots,n$.

All of the algorithms presented in this paper completely disentangle the structure of Problem (1) so that each iteration only involves the individual proximal operators of each of the nondifferentiable terms, the gradient operators of the differentiable terms, and multiplication by the linear maps. Thus, the maps $B_{i}$ are never inverted, and we never compute proximal operators or gradients of sums or infimal convolutions of functions. We note that this level of separability is not achieved by classical splitting methods such as forward-backward splitting, Douglas-Rachford splitting, or the alternating direction method of multipliers (ADMM) when they are applied directly to the primal optimization Problem (1) .

In Problem (1), the maps $B_{i}$ can be used as “data matrices,” in which case $h_{i}$ and $l_{i}$ are data fitting terms and $f$ and $g$ enforce prior knowledge on the structure of the solution, such as sparsity, low rank, or smoothness. In other cases, the maps $h_{i}$ and $l_{i}$ may be regularizers that emphasize many competing structures. We now present an example.

Application: Constrained model fitting with group-structured regularizers. Fix $d,m\in{\mathbf{N}}\backslash\{0\}$. Suppose we are given a measurement $b\in{\mathbf{R}}^{d}$ and a dictionary $A\in{\mathbf{R}}^{d\times m}$. Our goal is to recover a highly structured signal $x=(x_{1},\cdots,x_{m})^{T}\in{\mathbf{R}}^{m}$ such that $Ax\approx b$. For example, in the hierarchical sparse coding problem (HSCP) , we arrange the columns of $A$ into a directed tree structure ${\mathcal{T}}$ and allow $x_{i}=0$ only if $x_{j}=0$ for all descendants $j$ in ${\mathcal{T}}$ of node $i$. Such a hierarchical representation is particularly useful for multi-scale data such as images and text documents. This type of regularization can be generalized to include arbitrary column groupings and complicated relationships between the elements of each group. Indeed, let $G$ be a set of (possibly overlapping) subsets of $\{1,\cdots,m\}$. For all $S\in G$ and $x\in{\mathbf{R}}^{m}$, let $B_{S}x=L_{S}(x_{i})_{i\in S}^{T}\in{\mathbf{R}}^{m_{S}}$ where $m_{S}\in{\mathbf{N}}\backslash\{0\}$ and $L_{S}:{\mathbf{R}}^{|S|}\rightarrow{\mathbf{R}}^{m_{S}}$ is a linear map. Let $C\subseteq{\mathbf{R}}^{m}$ be a closed convex set, and let $\iota_{C}:{\mathbf{R}}^{m}\rightarrow\{0,\infty\}$ be the convex indicator function of $C$. For all $S\in G$, let $h_{S}:{\mathbf{R}}^{m_{S}}\rightarrow(-\infty,\infty]$ be a closed, proper, and convex regularizer, and let $l_{S}=\iota_{\{0\}}$, which implies $h_{S}\square l_{S}=h_{S}$. Then one special case of Problem (1) is the group-structured regularized model fitting problem:

Finally, we note that the use of infimal convolutions in applications is not wide-spread, so we list a few instances where they may be useful: Infimal convolutions are used in image recovery [14, Section 5] to remove staircasing effects in the total variation model. The infimal convolution of the indicator functions of two closed convex sets is the indicator function of their Minkowski sum, which has applications in motion planning for robotics [32, Section 4.3.2]. In convex analysis, the Moreau envelope of a function arises as an infimal convolution with a multiple of the squared norm [2, Section 12.4]. More generally, the infimal convolution of $h_{i}$ and $l_{i}$ can be interpreted as a regularization or smoothing of $h_{i}$ by $l_{i}$ and vice versa [2, Section 18.3].

This work seeks to improve the theoretical understanding of the convergence rates of primal-dual splitting schemes. In this paper, we study primal-dual algorithms that are applications of standard operator splitting algorithms in product spaces consisting of primal and dual variables. Consequently, the convergence theory for these algorithms is well-developed, and they are known to converge (weakly) under mild conditions.

Although we understand when these algorithms converge, relatively little is known about their rate of convergence. For convex optimization algorithms, the ergodic convergence rate of the primal-dual gap has been analyzed in a few cases . However, even in cases where convergence rates are known, variable metrics and stepsizes, which can significantly improve practical performance of the algorithms , are not analyzed. In addition, we are not aware of any convergence rate analysis of the primal-dual gap for the nonergodic (or last) iterate generated by these algorithms. It is important to understand nonergodic convergence rates because the ergodic (or time-averaged) iterates can “average out” structural properties, such as sparsity and low rank, that are shared by the solution and the nonergodic iterate.

The convergence rate analysis of the ergodic primal-dual gap largely follows from subgradient inequalities and an application of Jensen’s inequality. In contrast, the techniques developed in this paper exploit the properties of the nonexpansive operators driving the algorithms to deduce the nonergodic convergence rate of the primal-dual gap. Thus, our techniques are quite different from those used in classical convergence rate analysis and parallel the analysis developed in .

We summarize our contributions and techniques as follows:

We describe a model monotone inclusion problem that generalizes many primal-dual formulations that appear in the literature. We provide a simple prototype algorithm to solve the model problem, and we deduce a fundamental inequality that bounds the primal-dual gap at each iteration of the algorithm. We then simplify the inequality in the special case of four splitting algorithms (Section 2).

We derive ergodic convergence rates of the variable metric forms of the relaxed proximal point algorithm (PPA), relaxed forward-backward splitting (FBS), and forward-backward-forward splitting as well as the fixed metric relaxed Peaceman-Rachford splitting (PRS) algorithm (Section 3). After some algebraic simplifications, our analysis essentially follows from an application of Jensen’s inequality.

We derive nonergodic convergence rates of relaxed PPA, relaxed FBS, and relaxed PRS (Section 4). All of our analysis follows by bounding the primal-dual gap function by a multiple of the fixed-point residual (FPR) of the nonexpansive mapping that drives the algorithm. Thus, we show that the size of the FPR can be used as a valid stopping criteria for these three algorithms.

We apply our results to deduce ergodic and nonergodic convergence rates for a large class of primal-dual algorithms that have appeared in the literature (Section 5).

Our analysis not only deduces the convergence rates of a large class of primal-dual algorithms found in the literature. It also serves as a resource for the analysis of future primal-dual algorithms that solve generalizations of Problem 1, e.g., .

2 Definitions, notation and some facts

In what follows, ${\mathcal{H}},{\mathcal{G}}$, and ${\mathbf{H}}$ denote (possibly infinite dimensional) Hilbert spaces. We always use the notations $\langle\cdot,\cdot\rangle$ and $\|\cdot\|$ to denote the inner product and norm associated to a Hilbert space, respectively. Note that there is some ambiguity in this convention, but it simplifies the notation and no confusion should arise. The space ${\mathbf{H}}$ will usually denote a product Hilbert space consisting of primal variables in ${\mathcal{H}}$ and dual variables in ${\mathcal{G}}$. Let ${\mathbf{R}}_{++}=\{x\in{\mathbf{R}}\mid x>0\}$ denote the set of strictly positive real numbers. Let ${\mathbf{N}}=\{k\in{\mathbf{Z}}\mid k\geq 0\}$ denote the set of nonnegative integers. In all of the algorithms we consider, we utilize two stepsize sequences: the implicit sequence $(\gamma_{j})_{j\in{\mathbf{N}}}\subseteq{\mathbf{R}}_{++}$ and the explicit sequence $(\lambda_{j})_{j\in{\mathbf{N}}}\subseteq{\mathbf{R}}_{++}$. We define the $k$-th partial sum of the sequence $(\gamma_{j}\lambda_{j})_{j\in{\mathbf{N}}}$ by the formula:

Given a sequence $(x^{j})_{j\in{\mathbf{N}}}\subset{\mathcal{H}}$ and $k\in{\mathbf{N}}$, we let $\overline{x}^{k}=({1}/{\Sigma_{k}})\sum_{i=0}^{k}\gamma_{i}\lambda_{i}x^{i}$ denote its $k$th average with respect to the sequence $(\gamma_{j}\lambda_{j})_{j\in{\mathbf{N}}}$. We call a convergence result ergodic if it is in terms of the sequence $(\overline{x}^{j})_{j\in{\mathbf{N}}}$, and nonergodic if it is in terms of $(x^{j})_{j\in{\mathbf{N}}}$.

The following definitions and facts are mostly standard and can be found in

We let ${\mathcal{B}}({\mathcal{H}},{\mathcal{G}})$ denote the set of bounded linear maps from ${\mathcal{H}}$ to ${\mathcal{G}}$, and set ${\mathcal{B}}({\mathcal{H}}):={\mathcal{B}}({\mathcal{H}},{\mathcal{H}})$. We will use the notation $I_{{\mathcal{H}}}\in{\mathcal{B}}({\mathcal{H}})$ to denote the identity map. Given a map $L\in{\mathcal{B}}({\mathcal{H}},{\mathcal{G}})$, we denote its adjoint by $L^{\ast}\in{\mathcal{B}}({\mathcal{G}},{\mathcal{H}})$. The operator norm on $L\in{\mathcal{B}}({\mathcal{H}},{\mathcal{G}})$ is defined by the following supremum: $\|L\|=\sup_{x\in{\mathcal{H}},\|x\|\leq 1}\|Lx\|$. Let $\rho\in{\mathbf{R}}_{+}$ be a nonnegative real number. We let ${\mathcal{S}}_{\rho}({\mathcal{H}})\subseteq{\mathcal{B}}({\mathcal{H}})$ denote the set of linear $\rho$-strongly monotone self-adjoint maps:

We define the (semi)-norm and inner product induced by $U\in{\mathcal{S}}_{\rho}({\mathcal{H}})$ on ${\mathcal{H}}$ by the formulae: for all $x,y\in{\mathcal{H}}$, $\|x\|_{U}^{2}:=\langle Ux,x\rangle$, and $\langle x,y\rangle_{U}:=\langle Ux,y\rangle$. The Loewner partial ordering on ${\mathcal{S}}_{\rho}({\mathcal{H}})$ is defined as follows: for all $U_{1},U_{2}\in{\mathcal{S}}_{\rho}({\mathcal{H}})$, we have

Let $L\geq 0$, and let $D$ be a nonempty subset of ${\mathcal{H}}$. A map $T:D\rightarrow{\mathcal{H}}$ is called $L$-Lipschitz if for all $x,y\in D$, we have $\|Tx-Ty\|\leq L\|x-y\|$. In particular, $T$ is called nonexpansive if it is $1$-Lipschitz. A map $N:D\rightarrow{\mathcal{H}}$ is called $\lambda$-averaged [2, Section 4.4] if there exists a nonexpansive map $T:D\rightarrow{\mathcal{H}}$ and $\lambda\in(0,1)$ such that

A $(1/2)$-averaged map is called firmly nonexpansive.

Let $2^{\mathcal{H}}$ denote the power set of ${\mathcal{H}}$. A set-valued operator $A:{\mathcal{H}}\rightarrow 2^{\mathcal{H}}$ is called monotone if for all $x,y\in{\mathcal{H}}$, $u\in Ax$, and $v\in Ay$, we have $\langle x-y,u-v\rangle\geq 0$. We denote the set of zeros of a monotone operator by $\operatorname*{zer}(A):=\{x\in{\mathcal{H}}\mid 0\in Ax\}.$ The graph of $A$ is denoted by $\operatorname*{gra}(A):=\{(x,y)\mid x\in{\mathcal{H}},y\in Ax\}$. Evidently, $A$ is uniquely determined by its graph. A monotone operator $A$ is called maximal monotone provided that $\operatorname*{gra}(A)$ is not properly contained in the graph of any other monotone set-valued operator. The inverse of $A$, denoted by $A^{-1}$, is defined uniquely by its graph: $\operatorname*{gra}(A^{-1}):=\{(y,x)\mid x\in{\mathcal{H}},y\in Ax\}$. Let $\beta\in{\mathbf{R}}_{++}$ be a positive real number. The operator $A$ is called $\beta$-strongly monotone provided that for all $x,y\in{\mathcal{H}}$, $u\in Ax$, and $v\in Ay$, we have $\langle x-y,u-v\rangle\geq\beta\|x-y\|^{2}$. A single-valued operator $B:{\mathcal{H}}\rightarrow 2^{\mathcal{H}}$ maps each point in ${\mathcal{H}}$ to a singleton and will be identified with the natural ${\mathcal{H}}$-valued map it defines. A single-valued operator $B$ is called $\beta$-cocoercive provided that for all $x,y\in{\mathcal{H}}$, we have $\langle x-y,Bx-By\rangle\geq\beta\|Bx-By\|^{2}$. Evidently, $B$ is $\beta$-cocoercive whenever $B^{-1}$ is $\beta$-strongly monotone. The parallel sum of (not necessarily single-valued) monotone operators $A$ and $B$ is given by $A\square B:=(A^{-1}+B^{-1})^{-1}$. The resolvent of a monotone operator $A$ is defined by the inversion $J_{A}:=(I+A)^{-1}$. Minty’s theorem shows that $J_{A}$ is single-valued and has full domain ${\mathcal{H}}$ if, and only if, $A$ is maximally monotone. Note that $A$ is monotone if, and only if, $J_{A}$ is firmly nonexpansive. Thus, the reflection operator

is nonexpansive on ${\mathcal{H}}$ whenever $A$ is maximally monotone. If $\rho>0$ and $U\in{\mathcal{S}}_{\rho}({\mathcal{H}})$, the operator $U^{-1}A$ is maximal monotone in $\langle\cdot,\cdot\rangle_{U}$, if, and only if, $A$ is maximally monotone in $\langle\cdot,\cdot\rangle$. Let $\gamma\in(0,\infty)$. The resolvent of the map $\gamma U^{-1}A$ has the special identity: $J_{\gamma U^{-1}A}=U^{-1/2}J_{\gamma U^{-1/2}AU^{-1/2}}U^{1/2}$ [21, Example 3.9].

denote a subgradient of $f$ drawn at the point $x$, and the actual choice of the subgradient $\widetilde{\nabla}f(x)$ will always be clear from the context; note that this notation was also used in . The subdifferential operator of $f$ is maximally monotone. The inverse of $\partial f$ is given by $\partial f^{\ast}$ where $f^{\ast}(y):=\sup_{x\in{\mathcal{H}}}\{\langle y,x\rangle-f(x)\}$ is the Fenchel conjugate of $f$. If the function $f$ is $\beta$-strongly convex, then $\partial f$ is $\beta$-strongly monotone.

If a convex function $f:{\mathcal{H}}\rightarrow(-\infty,\infty]$ is Fréchet differentiable at $x\in{\mathcal{H}}$, then $\partial f(x)=\{\nabla f(x)\}$. Suppose $f$ is convex and Fréchet differentiable on ${\mathcal{H}}$, and let $\beta\in{\mathbf{R}}_{++}$ be a positive real number. Then the Baillon-Haddad theorem states that $\nabla f$ is $(1/\beta)$-Lipschitz, if, and only if, $\nabla f$ is $\beta$-cocoercive.

The resolvent operator associated to $\partial f$ is called the proximal operator and is uniquely defined by the following (strongly convex) minimization problem: $\mathbf{prox}_{f}(x):=J_{\partial f}(x)=\operatorname*{arg\,min}_{y\in{\mathcal{H}}}\{f(y)+(1/2)\|y-x\|^{2}\}$. If $\rho>0$, $U\in{\mathcal{S}}_{\rho}({\mathcal{H}})$, and $\gamma\in(0,\infty)$, the proximal operator of $f$ in the metric induced by $U$ is given by the following formula: for all $x\in{\mathcal{H}}$,

The infimal convolution of two functions $f,g:{\mathcal{H}}\rightarrow(-\infty,\infty]$ is denoted by $f\square g:{\mathcal{H}}\rightarrow[-\infty,\infty]:x\mapsto\inf_{y\in{\mathcal{H}}}\{f(y)+g(x-y)\}$. The indicator function of a closed, convex set $C\subseteq{\mathcal{H}}$ is denoted by $\iota_{C}:{\mathcal{H}}\rightarrow\{0,\infty\}$; the indicator function is on $C$ and is $\infty$ on ${\mathcal{H}}\backslash C$.

Finally, we call the following identity the cosine rule:

3 Assumptions

Every function we consider is closed, proper, and convex.

Unless otherwise stated, a function is not necessarily differentiable.

Every differentiable function we consider is Fréchet differentiable [2, Definition 2.45].

We employ other assumptions throughout the paper, but we list them closer to where they are invoked.

4 Basic properties of metrics

A simple proof of the following Lemma recently appeared in [20, Lemma 2.1]. It previously appeared in [30, Section VI.2.6].

Whenever $U,V\in{\mathcal{S}}_{0}({\mathcal{H}})$ satisfy the inequality $\alpha I_{{\mathcal{H}}}\succcurlyeq U\succcurlyeq V\succcurlyeq\beta I_{{\mathcal{H}}}$ for $\alpha,\beta>0$, we have the ordering $(1/\beta)I_{{\mathcal{H}}}\succcurlyeq V^{-1}\succcurlyeq U^{-1}\succcurlyeq(1/\alpha)I_{{\mathcal{H}}}$, the inclusion $U^{-1}\in{\mathcal{S}}_{\|U\|^{-1}}({\mathcal{H}})$, and the inequality $\|U^{-1}\|\leq(1/\beta)$.

5 Basic properties of resolvents and averaged operators

The following are simple modifications of standard facts found in .

Let $\rho>0$, let $\lambda>0$, let $\alpha\in(0,1)$, let $U\in{\mathcal{S}}_{\rho}({\mathcal{H}})$, let $A:{\mathcal{H}}\rightarrow{\mathcal{H}}$ be a single-valued maximal monotone operator, and let $f\in\Gamma_{0}({\mathcal{H}})$

Optimality conditions of $J$: We have $x^{+}:=J_{\gamma U^{-1}(\partial f+A)}(x)$ if, and only if, there exists a unique subgradient $\widetilde{\nabla}f(x^{+}):=(1/\gamma)U(x-x^{+})-Ax^{+}\in\partial f(x^{+})$, such that

Averaged operator contraction property: Let $\lambda\in(0,1)$. A map $T:{\mathcal{H}}\rightarrow{\mathcal{H}}$ is $\lambda$-averaged in the metric induced by $U$ if, and only if, for all $x,y\in{\mathcal{H}}$,

Wider relaxations: A map $T:{\mathcal{H}}\rightarrow{\mathcal{H}}$ is $\alpha$-averaged in $\|\cdot\|_{U}$, if, and only if, $T_{\lambda}$ (Equation (3)) is $\lambda\alpha$-averaged in $\|\cdot\|_{U}$ for all $\lambda\in(0,1/\alpha)$. In addition, $T_{1/\alpha}$ is nonexpansive with respect to $\|\cdot\|_{U}$.

6 Variable metrics

Throughout this paper we will consider sequences of mappings $(U_{j})_{j\in{\mathbf{N}}}\in{\mathcal{S}}_{\rho}({\mathcal{H}})$ for some $\rho>0$. In order to apply the standard convergence theory for variable metrics, we will make the following assumption:

Assumption 3 is standard in variable metric algorithms .

There is an asymmetry in our notation and the notation of . In our analysis, the map $U\in{\mathcal{S}}_{\rho}({\mathcal{H}})$ induces a metric on ${\mathcal{H}}$. In other papers, the maps $U^{-1}$ induce a metric on ${\mathcal{H}}$.

The following notation will be used throughout the rest of the paper. The proof is elementary.

The following Proposition is a consequence of the proof of [20, Theorem 5.1]. The proof is simple, so we omit it.

For all $k\in{\mathbf{N}}$, $\|z^{k+1}-z^{\ast}\|_{U_{k+1}}^{2}\leq(1+\eta_{k})\|z^{k}-z^{\ast}\|_{U_{k}}^{2}$ and hence,

We will use the following proposition to select parameters in the FBS algorithm. The proof of the following fact follows from [46, Equation (3.35)]

Let $\rho>0$, let $\beta>0$, let $B:{\mathcal{H}}\rightarrow{\mathcal{H}}$ be $\beta$-cocoercive in the norm $\|\cdot\|$, and let $U\in{\mathcal{S}}_{\rho}({\mathcal{H}})$. Then $U^{-1}B$ is $\beta\rho$-cocoercive in the norm $\|\cdot\|_{U}$.

The following proposition essentially follows from the proof of [45, Theorem 3.1].

Let $A:{\mathcal{H}}\rightarrow 2^{\mathcal{H}}$ be maximal monotone, let $B:{\mathcal{H}}\rightarrow{\mathcal{H}}$ be monotone and $(1/\beta)$-Lipschitz for some $\beta>0$, let $\rho>0$, let $(U_{j})_{j\in{\mathbf{N}}}\subseteq{\mathcal{S}}_{\rho}({\mathcal{H}})$ satisfy Assumption 3, and let $(\gamma_{j})_{j\in{\mathbf{N}}}\subseteq(0,\rho\beta]$. Let $(z^{j})_{j\in{\mathbf{N}}}$ be a sequence of points defined by the iteration: let $z^{0}\in{\mathcal{H}}$ and for all $k\in{\mathbf{N}}$, define

Suppose that $\operatorname*{zer}(A+B)\neq\emptyset$. Then for all $z^{\ast}\in\operatorname*{zer}(A+B)$ and for all $k\in{\mathbf{N}}$, we have, $\|z^{k+1}-z^{\ast}\|_{U_{k+1}}^{2}\leq(1+\eta_{k})\|z^{k}-z^{\ast}\|_{U_{k}}^{2}$.

The unifying scheme

In this section, we introduce a prototype monotone inclusion problem that generalizes and summarizes many primal-dual problem formulations found in the literature. After we describe the problem, we will introduce an abstract unifying scheme that generalizes many existing primal-dual algorithms. We will describe how to measure convergence of the unifying scheme, and introduce a fundamental inequality that bounds our measure of convergence. Finally, we will identify the key terms in the fundamental inequality and simplify them in the case of several abstract splitting algorithms.

In Section 5, we will show that this unifying scheme relates to many existing algorithms, and extend the convergence rate results of those methods.

Let $({\mathbf{H}},\langle\cdot,\cdot\rangle)$ be a Hilbert space, let ${\mathbf{f}},{\mathbf{g}}\in\Gamma_{0}({\mathbf{H}})$, and let ${\mathbf{S}}:{\mathbf{H}}\rightarrow{\mathbf{H}}$ be a skew symmetric map: ${\mathbf{S}}^{\ast}=-{\mathbf{S}}$. Then the prototype primal-dual problem is to find ${\mathbf{x}}^{\ast}\in{\mathbf{H}}$ such that

Evidently, Problem 1 is a monotone inclusion problem because $\partial{\mathbf{f}},\partial{\mathbf{g}}$, and ${\mathbf{S}}$ are maximally monotone operators on ${\mathbf{H}}$ [2, Example 20.30].

We are now ready to define our unifying scheme.

Note that the points ${\mathbf{x}}_{\mathbf{f}}^{k},{\mathbf{x}}_{\mathbf{g}}^{k}$, and ${\mathbf{x}}_{\mathbf{S}}^{k}$ as well as the subgradients $\widetilde{\nabla}{\mathbf{f}}({\mathbf{x}}_{\mathbf{f}}^{k})\in\partial{\mathbf{f}}({\mathbf{x}}_{\mathbf{f}}^{k})$ and $\widetilde{\nabla}{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}}^{k})\in\partial{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}}^{k})$ are unspecified in the description of Algorithm 1. In the algorithms we study, these points and subgradients will be generated by proximal and forward gradient operators and, thus, can be determined given ${\mathbf{z}}^{k}$; see Section 2.2 for examples. However, Algorithm 1 is only meant to illustrate the algebraic form that our analysis addresses, and it is not meant to be an actual algorithm that solves Problem 9. The positive scalar sequence $(\lambda_{j})_{j\in{\mathbf{N}}}$ consists of relaxation parameters, or explicit stepsize parameters, whereas the sequence $(\gamma_{j})_{j\in{\mathbf{N}}}$ consists of proximal parameters, or implicit stepsize parameters. The strongly monotone maps $(U_{j})_{j\in{\mathbf{N}}}$ induce the metrics used in each iteration of the algorithm.

In all of our applications, ${\mathbf{H}}$ will be a product space of primal and dual variables. In this setting, ${\mathbf{f}}$ and ${\mathbf{g}}$ will be block-separable maps, and ${\mathbf{g}}$ will sometimes be differentiable. The map ${\mathbf{S}}$ “mixes” the primal and dual variable sequences in the product space. Mixing is necessary, because the sequences are otherwise uncoupled.

The sequence of maps $(U_{j})_{j\in{\mathbf{N}}}$ is employed for two purposes. First, the maps are used because the evaluation of the resolvent $J_{\partial{\mathbf{f}}+{\mathbf{S}}}$, which is a basic building block of most of the algorithms we study, may not be simple. Thus, the primal-dual algorithms that we study formulate special metrics induced by $U\in{\mathcal{S}}_{\rho}({\mathbf{H}})$ such that $J_{U^{-1}(\partial{\mathbf{f}}+{\mathbf{S}})}$ is as easy to evaluate as $\mathbf{prox}_{{\mathbf{f}}}$ (See Section 5). Hence, in our analysis we must at least consider fixed metrics that are different from the standard product metric on ${\mathbf{H}}$. Second, we allow the metrics to vary at each iteration because it can significantly improve the practical performance of the algorithm, e.g., by employing second order information, or even simple time-varying diagonal metrics .

2 Examples of the unifying scheme

In this section we introduce four algorithms and show that they are special cases of Algorithm 1. We will also introduce several assumptions on the algorithm parameters that ensure convergence. These assumptions will remain in effect throughout the rest of the paper. Note that the convergence theory of the methods in this section is well-studied. See for background. Finally, we will say that several algorithms in this section are relaxed. For brevity, we will drop this adjective whenever convenient.

The relaxed variable metric PPA applies to problems in which ${\mathbf{g}}\equiv 0$.

The relaxed variable metric FBS algorithm can be applied whenever ${\mathbf{g}}$ is differentiable and $\nabla{\mathbf{g}}$ is $({1}/{\beta})$-Lipschitz for some $\beta>0$.

In the relaxed PRS algorithm, we fix the metric and the implicit stepsize parameters throughout the course of the algorithm. We do this because the fixed-points of the PRS operator can vary with $\gamma$ and $U$. Thus, changing these parameters will lead to an algorithm that “chases” a new fixed-point at each iteration.

The variable metric FBF algorithm can be applied whenever ${\mathbf{g}}$ is differentiable and $\nabla{\mathbf{g}}$ is $({1}/{\beta})$-Lipschitz for some $\beta>0$.

The following lemma relates the above algorithms to the unifying scheme.

Algorithms 2, 3, 4, and 5 are special cases of the unifying scheme. In particular, the following hold for all $k\in{\mathbf{N}}$:

In Algorithm 2, we have ${\mathbf{x}}_{\mathbf{g}}^{k}:={\mathbf{z}}^{k}$, ${\mathbf{x}}_{\mathbf{S}}^{k}:={\mathbf{x}}_{\mathbf{f}}^{k},$ and

In Algorithm 3, we have ${\mathbf{x}}_{\mathbf{g}}^{k}:={\mathbf{z}}^{k}$, ${\mathbf{x}}_{\mathbf{S}}^{k}:={\mathbf{x}}_{\mathbf{f}}^{k}$, and

In Algorithm 4, we have ${\mathbf{z}}^{k+1}-{\mathbf{z}}^{k}=\lambda_{k}({\mathbf{x}}_{{\mathbf{f}}}^{k}-{\mathbf{x}}_{{\mathbf{g}}}^{k})$ for

and $\widetilde{\nabla}{\mathbf{f}}({\mathbf{x}}_{\mathbf{f}}^{k}):=(1/\gamma)U(2{\mathbf{x}}_{\mathbf{g}}^{k}-{\mathbf{z}}^{k}-{\mathbf{x}}_{\mathbf{f}}^{k})-w{\mathbf{S}}{\mathbf{x}}_{\mathbf{f}}^{k}\in\partial{\mathbf{f}}({\mathbf{x}}_{\mathbf{f}}^{k})$.

In Algorithm 5, we have $\lambda_{k}=1$, ${\mathbf{x}}_{\mathbf{g}}^{k}:={\mathbf{x}}_{\mathbf{f}}^{k}$, ${\mathbf{x}}_{\mathbf{S}}^{k}:={\mathbf{x}}_{\mathbf{f}}^{k}$, and

Proof. Fix $k\in{\mathbf{N}}$, and note that the subgradient identities all follow from Part 1 of Proposition 2.

Part 2: From Part 1 of Proposition 2, we have the following identity:

Thus, altogether we have ${\mathbf{z}}^{k+1}={\mathbf{z}}^{k}-\gamma_{k}\lambda_{k}U_{k}^{-1}\left(\widetilde{\nabla}{\mathbf{f}}({\mathbf{x}}_{\mathbf{f}}^{k})+\nabla{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}}^{k})+{\mathbf{S}}{\mathbf{x}}_{{\mathbf{S}}}^{k}\right).$

Therefore, if we define ${\mathbf{x}}_{\mathbf{S}}^{k}:=w{\mathbf{x}}_{\mathbf{f}}^{k}+(1-w){\mathbf{x}}_{\mathbf{g}}^{k}$, then

Now we establish two basic and well known results on the boundedness and summability of various terms related to the above algorithms. These facts will be used repeatedly in our convergence rate analysis.

Let $\rho\in{\mathbf{R}}_{++}$, let $U\in{\mathcal{S}}_{\rho}({\mathbf{H}})$, let $\gamma\in{\mathbf{R}}_{++}$, and let $\beta\in{\mathbf{R}}_{++}$. Then the following hold:

The operator $J_{\gamma U^{-1}(\partial{\mathbf{f}}+{\mathbf{S}})}$ is $(1/2)$-averaged in the norm $\|\cdot\|_{U}$. In addition, the set of fixed points of $J_{\gamma U^{-1}(\partial{\mathbf{f}}+{\mathbf{S}})}$ is equal to $\operatorname*{zer}\left(\partial{\mathbf{f}}+{\mathbf{S}}\right)$

Let $\gamma\in(0,2\beta\rho)$. Suppose that ${\mathbf{g}}$ is differentiable and $\nabla{\mathbf{g}}$ is $(1/\beta)$-Lipschitz. Then the composition

is $\alpha_{\rho,\gamma}$-averaged in the norm $\|\cdot\|_{U}$ where

Let $w\in{\mathbf{R}}$, and define the PRS operator:

Parts 1 and 3 are simple modifications of standard facts found in .

Part 2: Note that $U^{-1}\nabla{\mathbf{g}}$ is $\beta\rho$-cocoercive in $\|\cdot\|_{U}$ by Proposition 5 and the Baillon-Haddad theorem . Thus, $I_{{\mathbf{H}}}-\gamma U^{-1}\nabla{\mathbf{g}}$ is $\gamma/(2\beta\rho)$ averaged in $\|\cdot\|_{U}$ by [2, Proposition 4.33]. Thus, the formula for $\alpha_{\rho,\gamma}$ follows from [36, Theorem 3(b)]. The fixed-point identity follows from a simple modification of [2, Theorem 25.1]. ∎

Let ${\mathbf{z}}^{\ast}\in\operatorname*{zer}(\partial{\mathbf{f}}+\nabla{\mathbf{g}}+{\mathbf{S}})$. Then in Algorithm 3, the following are true:

Part 4 follows from Proposition 6 applied to the the maximal monotone operator $\partial{\mathbf{f}}$ and the $(\beta^{-1}+\|{\mathbf{S}}\|)$-Lipschitz operator $\nabla{\mathbf{g}}+{\mathbf{S}}$. ∎

3 The fundamental inequality

This section describes the pre-primal-dual gap (Definition 11). We use the pre-primal-dual gap to measure the convergence of the unifying scheme. In Section 5, we will show that under certain conditions, the pre-primal-dual gap function bounds the primal and dual objective errors of the iterates generated by a class of primal-dual algorithms.

Before we introduce the gap function, we analyze the optimality conditions of Problem 1. The following lemma is well-known.

Let ${\mathbf{x}}^{\ast}\in{\mathbf{H}}$. Suppose that ${\mathbf{x}}^{\ast}$ solves Problem 1. Then for all ${\mathbf{x}}\in{\mathbf{H}},$

If ${\mathbf{x}}^{\ast}$ solves Problem 1, then $-{\mathbf{S}}{\mathbf{x}}^{\ast}$ is a subgradient of ${\mathbf{f}}+{\mathbf{g}}$ at the point ${\mathbf{x}}^{\ast}$. Thus, Equation (15) follows after noting that $\langle{\mathbf{S}}{\mathbf{x}},{\mathbf{x}}\rangle=0$ for all ${\mathbf{x}}\in{\mathbf{H}}$.

The other direction follows because Equation (15) characterizes the set of subgradients of the form $-{\mathbf{S}}{\mathbf{x}}^{\ast}\in\partial({\mathbf{f}}+{\mathbf{g}})({\mathbf{x}}^{\ast})=\partial{\mathbf{f}}({\mathbf{x}}^{\ast})+\partial{\mathbf{g}}({\mathbf{x}}^{\ast})$. ∎

See [2, Corollary 16.38] for conditions that imply additivity of the subdifferential.

Lemma 10 motivates the following definition:

Let the setting be as in Algorithm 1. Define the pre-primal dual gap function by the formula: for all ${\mathbf{x}}_{\mathbf{f}},{\mathbf{x}}_{\mathbf{g}},{\mathbf{x}}_{{\mathbf{S}}},{\mathbf{x}}\in{\mathbf{H}}$, let

then ${\mathbf{x}}^{\prime}$ is a solution of Problem 1 (Lemma 10).

Finally, Lemma 10 shows that for all ${\mathbf{x}}\in{\mathbf{H}}$,

whenever ${\mathbf{x}}^{\ast}$ solves Problem 1. See Section 5.1 for other lower bounds of the pre-primal-dual gap in the context of a particular convex optimization problem.

The following is our main tool to bound the pre-primal-dual gap.

Suppose that $({\mathbf{z}}^{j})_{j\geq 0}$ is generated by Algorithm 1, and let ${\mathbf{x}}\in{\mathbf{H}}$. Then the following inequality holds: for all $k\in{\mathbf{N}}$,

Fix $k\in{\mathbf{N}}$. First expand the norm:

We add and subtract a point in the inner products involving ${\mathbf{f}}$ and ${\mathbf{g}}$ and use the subgradient inequality to get:

Therefore Equation (19) follows after rearranging. ∎

The upper fundamental inequality in Proposition 12 bounds the pre-primal-dual gap with the sum of an alternating sequence and a key term.

Let $({\mathbf{z}}^{j})_{j\in{\mathbf{N}}}$ be generated by Algorithm 1. For all $k\in{\mathbf{N}}$, we define the fundamental upper key term

The value $\kappa_{u}^{k}(\lambda_{k})$ depends on the entire history of Algorithm 1 up to and including iteration $k$, but in our analysis we will only view $\kappa_{u}^{k}(\lambda_{k})$ as a function of the parameter $\lambda_{k}$. Throughout the rest of the paper, we will often make the dependence of the upper key term on $\lambda_{k}$ implicit, and denote $\kappa_{u}^{k}:=\kappa_{u}^{k}(\lambda_{k})$. However, in the proof of Theorem 4 we will need to keep the dependence explicit.

The following proposition will compute the upper key terms induced by the PPA, FBS, PRS, and FBF algorithms. See Section 2.2 for the definitions of the points ${\mathbf{x}}_{\mathbf{f}}^{k},{\mathbf{x}}_{\mathbf{g}}^{k}$, and ${\mathbf{x}}_{\mathbf{S}}^{k}$.

Let $({\mathbf{z}}^{j})_{j\in{\mathbf{N}}}$ be generated by Algorithm 1. Then for all $k\in{\mathbf{N}}$, the following inequalities and identities hold:

In Algorithm 2, we have $\kappa_{u}^{k}(\lambda_{k})=\left(1-2/\lambda_{k}\right)\|{\mathbf{z}}^{k+1}-{\mathbf{z}}^{k}\|_{U_{k}}^{2}.$

In Algorithm 4, we have $\kappa_{u}^{k}(\lambda_{k})=\left(1-2/\lambda_{k}\right)\|{\mathbf{z}}^{k+1}-{\mathbf{z}}^{k}\|_{U}^{2}.$

In Algorithm 5, we have $\kappa_{u}^{k}(\lambda_{k})\leq 0.$

Proof. Fix $k\in{\mathbf{N}}$. To simplify notation, we drop the iteration index and denote ${\mathbf{z}}:={\mathbf{z}}^{k},{\mathbf{x}}_{\mathbf{f}}:={\mathbf{x}}_{\mathbf{f}}^{k},{\mathbf{x}}_{\mathbf{g}}:={\mathbf{x}}_{\mathbf{g}}^{k},{\mathbf{x}}_{\mathbf{S}}:={\mathbf{x}}_{\mathbf{S}}^{k},{\mathbf{z}}^{+}:={\mathbf{z}}^{k+1},\gamma:=\gamma_{k},\lambda:=\lambda_{k},U:=U_{k},$ and $\kappa_{u}:=\kappa_{u}^{k}(\lambda_{k})$ throughout this proof.

For PPA, FBS, and PRS, we note that the following identities hold:

and there exists $w\in{\mathbf{R}}$ such that

Indeed, in PPA and FBS, $w=1$ (see Section 2.2). In PRS, $w$ is a parameter of the algorithm, and Equations (22) and (21) are shown in Lemma 7. Furthermore, Part 1 of Proposition 2 shows that in PPA and FBS,

for a unique subgradient $\widetilde{\nabla}{\mathbf{f}}({\mathbf{x}}_{\mathbf{f}})\in\partial{\mathbf{f}}({\mathbf{x}}_{\mathbf{f}})$; see Lemma 7 for the definition of $\widetilde{\nabla}{\mathbf{f}}({\mathbf{x}}_{\mathbf{f}})$.

where we make the identification $\widetilde{\nabla}{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}})=\nabla{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}})$ whenever ${\mathbf{g}}$ is differentiable; see Lemma 7 for the definition of $\widetilde{\nabla}{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}})$ in the PRS algorithm. Because ${\mathbf{x}}_{\mathbf{S}}={\mathbf{x}}_{\mathbf{g}}+w({\mathbf{x}}_{\mathbf{f}}-{\mathbf{x}}_{\mathbf{g}})={\mathbf{x}}_{\mathbf{g}}+(w/\lambda)({\mathbf{z}}^{+}-{\mathbf{z}})$ and $\langle{\mathbf{S}}{\mathbf{x}},{\mathbf{x}}\rangle=0$ for all ${\mathbf{x}}\in{\mathbf{H}}$, we have the simplification:

where the second to last equality uses Equation (23) and the second to last “$+$” also uses Equation (22).

Now we proceed with the specific cases: In PPA and FBS, $w=1$ and

where use the identity ${\mathbf{x}}_{\mathbf{f}}-{\mathbf{z}}^{+}=\left(1-(1/\lambda)\right)({\mathbf{z}}-{\mathbf{z}}^{+})$ on the third line, we use the identity ${\mathbf{z}}^{+}-{\mathbf{z}}=\lambda({\mathbf{x}}_{\mathbf{f}}-{\mathbf{x}}_{\mathbf{g}})$ (Equation (21)) on the last two lines, and the last inequality follows from the Descent Theorem [2, Theorem 18.15(iii)]: $\langle\nabla{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}}),{\mathbf{x}}_{\mathbf{g}}-{\mathbf{x}}_{\mathbf{f}}\rangle\leq{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}})-{\mathbf{g}}({\mathbf{x}}_{\mathbf{f}})+(1/(2\beta))\|{\mathbf{x}}_{\mathbf{g}}-{\mathbf{x}}_{\mathbf{f}}\|^{2}.$ In PPA ${\mathbf{g}}\equiv 0$, so the Equation (26) implies the identity in Part 1. The inequality for FBS now follows by the above bound for $\kappa_{u}$, the bound $\gamma\leq 2\beta\rho-\varepsilon$, and

where we use $\lambda\leq(4\beta\rho-\gamma)/2\beta\rho\leq 2$ and the lower bound $U\succcurlyeq\rho I_{{\mathbf{H}}}$.

Therefore, subtract $\lambda{\mathbf{z}}^{+}+(1-\lambda){\mathbf{z}}$ from both sides of the above equation, divide by $\lambda$, and use the identity in Equation (24) to get

Finally, we prove the bound for the FBF algorithm:

Note that the operator $\nabla{\mathbf{g}}+{\mathbf{S}}$ is $(1/\beta)+\|{\mathbf{S}}\|$ Lipschitz. Thus,

where we use the following bound: for all ${\mathbf{x}}\in{\mathbf{H}}$, $\|{\mathbf{x}}\|_{U^{-1}}^{2}\leq(1/\rho)\|{\mathbf{x}}\|^{2}$ (Lemma 1).

Ergodic convergence

In this section, we prove an ergodic convergence rate for the pre-primal-dual gap. To this end, we recall the partial sum sequence $\Sigma_{k}=\sum_{i=0}^{k}\gamma_{i}\lambda_{i},$ and for every sequence of vectors $({\mathbf{x}}^{j})_{j\geq 0}\subseteq{\mathbf{H}}$, we define the ergodic sequence $\overline{{\mathbf{x}}}^{k}=({1}/{\Sigma_{k}})\sum_{i=0}^{k}\gamma_{i}\lambda_{i}{\mathbf{x}}^{i}$. For each algorithm, Theorem 16 (below) gives an ergodic sequence $(\overline{{\mathbf{x}}}^{j})_{j\in{\mathbf{N}}}$ such that for all bounded subsets $D\subseteq{\mathbf{H}}$, we have

This bound is a generalization of the primal-dual gap bounds shown in . See Section 5.1 for several lower bounds of the pre-primal-dual gap.

Before we prove our ergodic rates, we need to prove a bound for PRS. Recall that we only analyze the PRS algorithm when the map $U_{k}\equiv U$ is fixed. The following lemma will help us deduce the convergence rate of the PRS algorithm whenever ${\mathbf{f}}$ or ${\mathbf{g}}$ is Lipschitz (Part 3 of Theorem 16).

Proof. Fix $k\in{\mathbf{N}}$. The identity $\lambda_{k}({\mathbf{x}}_{\mathbf{f}}^{k}-{\mathbf{x}}_{{\mathbf{g}}}^{k})={\mathbf{z}}^{k+1}-{\mathbf{z}}^{k}$ and the fact the sequence $(\|{\mathbf{z}}^{j}-{\mathbf{z}}^{\ast}\|_{U})_{j\in{\mathbf{N}}}$ is decreasing (Part 3 of Proposition 9), show that

Lemma 15 shows that the difference of splitting variables $\overline{{\mathbf{x}}}_{{\mathbf{f}}}^{k}-\overline{{\mathbf{x}}}_{\mathbf{g}}^{k}$ converges to zero with rate $O({1}/{\Sigma_{k}})$. Thus, if ${\mathbf{f}}$ is Lipschitz continuous, then $|{\mathbf{f}}(\overline{{\mathbf{x}}}_{\mathbf{f}}^{k})-{\mathbf{f}}(\overline{{\mathbf{x}}}_{\mathbf{g}}^{k})|=O(1/\Sigma_{k})$.

We are now ready to prove our main ergodic convergence results.

Suppose that the sequence $({\mathbf{z}}^{j})_{j\in{\mathbf{N}}}$ is generated by Algorithm 1, and suppose that Assumption 3 holds. Then for all ${\mathbf{x}}\in{\mathbf{H}}$ and all $k\in{\mathbf{N}}$, we have the following bounds:

Ergodic convergence of PPA: Let ${\mathbf{z}}^{\ast}\in\operatorname*{zer}(\partial{\mathbf{f}}+{\mathbf{S}})$. Then in Algorithm 2, we have

Ergodic convergence of FBS: Let ${\mathbf{z}}^{\ast}\in\operatorname*{zer}(\partial{\mathbf{f}}+\nabla{\mathbf{g}}+{\mathbf{S}})$, and let $\overline{\lambda}=\sup_{j\in{\mathbf{N}}}\lambda_{j}$. Then in Algorithm 3, we have the bounds $0<\inf_{j\in{\mathbf{N}}}\lambda_{j}\leq\overline{\lambda}\leq 2$ and $\inf_{j\in{\mathbf{N}}}(1-\alpha_{j}\lambda_{j})/(\alpha_{j}\lambda_{j})>0$, and

Ergodic convergence of FBF: Let ${\mathbf{z}}^{\ast}\in\operatorname*{zer}(\partial{\mathbf{f}}+\nabla{\mathbf{g}}+{\mathbf{S}})$. Then in Algorithm 5, we have

Proof. Fix $k\in{\mathbf{N}}$. For any sequence of points $({\mathbf{z}}^{j})_{j\in{\mathbf{N}}}\subseteq{\mathbf{H}}$ and any point ${\mathbf{z}}^{\ast}\in{\mathbf{H}}$ such that $\|{\mathbf{z}}^{i+1}-{\mathbf{z}}^{\ast}\|_{U_{i+1}}^{2}\leq(1+\eta_{i})\|{\mathbf{z}}^{i}-{\mathbf{z}}^{\ast}\|_{U_{i}}^{2}$ for all $i\in{\mathbf{N}}$, we have $\|{\mathbf{z}}^{i}-{\mathbf{z}}^{\ast}\|^{2}_{U_{i}}\leq\left(\prod_{i=0}^{\infty}(1+\eta_{i})\right)\|{\mathbf{z}}^{0}-{\mathbf{z}}^{\ast}\|_{U_{0}}^{2}.$ Therefore, by the convexity of $\|\cdot\|_{U_{i}}^{2}$ for all $i\in{\mathbf{N}}$, and by the inequality $-\|{\mathbf{x}}\|_{U_{i}}\leq-(1/(1+\eta_{i}))\|{\mathbf{x}}\|_{U_{i+1}}$ for all ${\mathbf{x}}\in{\mathbf{H}}$ and $i\in{\mathbf{N}}$, we have

We will use Equation (29) to produce bounds for all of the variable metric methods.

Part 2: We have the following bound from Proposition 9:

Thus, the bound follows from Jensen’s inequality, Proposition 14 ($\kappa_{u}^{i}\leq\left(\rho-{\varepsilon}/{(\beta\lambda_{i})}\right)\|{\mathbf{z}}^{i+1}-{\mathbf{z}}^{i}\|^{2}+2\gamma_{i}\lambda_{i}{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}}^{i})-2\gamma_{i}\lambda_{i}{\mathbf{g}}({\mathbf{x}}_{\mathbf{f}}^{i})$), and the fundamental inequality:

Part 3: We prove the result when ${\mathbf{f}}$ is Lipschitz; the other case is symmetric. This follows from the Jensen’s inequality, Proposition 14 ($\kappa_{u}^{i}=\left(1-{2}/{\lambda_{i}}\right)\|{\mathbf{z}}^{i+1}-{\mathbf{z}}^{i}\|_{U}^{2}\leq 0$), the fundamental inequality, and the identity $\overline{{\mathbf{x}}}_{\mathbf{g}}^{k}-\overline{{\mathbf{x}}}_{\mathbf{S}}^{k}=w(\overline{{\mathbf{x}}}_{\mathbf{g}}^{k}-\overline{{\mathbf{x}}}_{\mathbf{f}}^{k})$ (follows by averaging identities found in Part 3 of Lemma 7):

Part 4: This follows from the Jensen’s inequality, Proposition 14 ($\kappa_{u}^{i}\leq 0$), and the fundamental inequality:

In general, the $O(1/(k+1))$ convergence rates in Theorem 16 are the best PPA, FBS, and PRS obtain for $(\overline{{\mathbf{x}}}_{\mathbf{f}}^{j})_{j\in{\mathbf{N}}}$ and $(\overline{{\mathbf{x}}}_{{\mathbf{g}}}^{j})_{j\in{\mathbf{N}}}$ [24, Proposition 8].

Nonergodic convergence

In this section we deduce nonergodic convergence rates for PPA, FBS and PRS under the following assumption:

For all nonergodic convergence results, we assume $(U_{j})_{j\in{\mathbf{N}}}$ and $(\gamma_{j})_{j\in{\mathbf{N}}}$ are constant sequences.

For PPA, FBS, and PRS, Theorem 18 (below) produces a natural sequence $({\mathbf{x}}^{j})_{j\in{\mathbf{N}}}$ such that for all bounded subsets $D\subseteq{\mathbf{H}}$, we have

To the best of our knowledge, the rate of convergence for the nonergodic primal-dual gap generated by the class of algorithms we study has never appeared in the literature.

Nonergodic iterates tend to share structural properties, such as sparsity or low rank, with the solution of the problem. In some cases, the ergodic iterates generated in Section 3 “average out” structural properties of the nonergodic iterates. Thus, although the ergodic iterates may be “closer” to the solution, they are often poorer partial solutions than the nonergodic iterates. The results of this section provide worst-case theoretical guarantees on the quality of the nonergodic iterates in order to justify their use in practical applications.

In our analysis, we use the following result (see also for similar little-$o$ and big-$O$ convergence rates):

Let $\alpha\in(0,1)$, let $\rho>0$, let $U\in{\mathcal{S}}_{\rho}({\mathbf{H}})$, and let $(\lambda_{j})_{j\in{\mathbf{N}}}\subseteq(0,1/\alpha)$. Suppose that $T:{\mathbf{H}}\rightarrow{\mathbf{H}}$ is an $\alpha$-averaged operator in the norm $\|\cdot\|_{U}$. Let ${\mathbf{z}}^{\ast}$ be a fixed point of $T$, let ${\mathbf{z}}^{0}\in{\mathbf{H}}$, let $\tau_{k}:=(1-\alpha\lambda_{k})\lambda_{k}/\alpha$ for all $k\in{\mathbf{N}}$, suppose that $\underline{\tau}:=\inf_{j\in{\mathbf{N}}}\tau_{j}>0$, and suppose that $({\mathbf{z}}^{j})_{j\in{\mathbf{N}}}$ is generated by the following iteration: for all $k\in{\mathbf{N}}$, let

Throughout this section, $T$ will always denote an $\alpha$-averaged mapping in the norm $\|\cdot\|_{U}$. Recall that for $\lambda\in(0,1/\alpha)$, $T_{\lambda}$ is $\alpha\lambda$-averaged (see Proposition 2), so

for all $k\in{\mathbf{N}}$, and any fixed-point ${\mathbf{z}}^{\ast}$ of $T$. Note that Equation (33) also holds when $\alpha\lambda=1$ (see Proposition 2). Equation (33) shows that $T_{\lambda}{\mathbf{z}}^{k}$ is at least as close to ${\mathbf{z}}^{\ast}$ as ${\mathbf{z}}^{k}$ is. This fact will be useful in the proof of Theorem 18 below.

In the following theorem, we will deduce little-$o$ and big-$O$ convergence rates. Because the pre-primal-dual gap can be negative, we slightly abuse notation: given a point ${\mathbf{x}}\in{\mathbf{H}}$, a (not necessarily positive) sequence $(a_{j})_{j\in{\mathbf{N}}}$ satisfies $a_{k}=o((1+\|{\mathbf{x}}\|_{U})/\sqrt{k+1})$ provided that there exists a nonnegative sequence $(b_{j})_{j\in{\mathbf{N}}}$ such that $b_{k}=o((1+\|{\mathbf{x}}\|_{U})/\sqrt{k+1})$ and $a_{k}=O(b_{k})$. Note that we do not measure $|a_{k}|$ because our only goal is to ensure that the sequence $(a_{j})_{j\in{\mathbf{N}}}$ is eventually nonpositive.

Suppose that Assumption 4 holds, let $U\in{\mathcal{S}}_{\rho}({\mathbf{H}})$ denote the common metric inducing map, and let $\gamma\in{\mathbf{R}}_{++}$ denote the common stepsize parameter. Then each method is a special case of Iteration (31). For each method, assume that $\underline{\tau}>0$ (See Theorem 17). Then for all $k\in{\mathbf{N}}$ and all ${\mathbf{x}}\in{\mathbf{H}}$, the following hold:

Nonergodic convergence of PPA: Let ${\mathbf{z}}^{\ast}\in\operatorname*{zer}(\partial{\mathbf{f}}+{\mathbf{S}})$. Then in Algorithm 2, we have $\alpha=1/2$ and $T=J_{U^{-1}(\partial{\mathbf{f}}+{\mathbf{S}})}$,

Proof. Fix $k\in{\mathbf{N}}$. In all of the following proofs, we will bound the pre-primal-dual gap by a quantity involving $\|T{\mathbf{z}}^{k}-{\mathbf{z}}^{k}\|_{U}$. Then the big-$O$ and little-$o$ convergence rates follow directly from Theorem 17. In addition, we will use Equation (33) and the independence of ${\mathbf{x}}_{\mathbf{f}}^{k},{\mathbf{x}}_{\mathbf{g}}^{k},$ and ${\mathbf{x}}_{\mathbf{S}}^{k}$ from $\lambda_{k}$ to tighten our upper bounds. To this end, we will denote ${\mathbf{z}}_{\lambda}:=T_{\lambda}({\mathbf{z}}^{k})$ (see Equation (3)) and let $C=(0,1/\alpha]$ where $\alpha$ is averagedness coefficient of $T$. Note that $T_{\lambda}$ is nonexpansive for all $\lambda\in C$ (see Part 3 of Proposition 2). Also note that for $\lambda\in C$, we have $(1/\lambda)({\mathbf{z}}_{\lambda}-{\mathbf{z}}^{k})=T{\mathbf{z}}^{k}-{\mathbf{z}}^{k}$ and $\|{\mathbf{z}}_{\lambda}-{\mathbf{z}}^{\ast}\|_{U}\leq\|{\mathbf{z}}^{k}-{\mathbf{z}}^{\ast}\|_{U}\leq\|{\mathbf{z}}^{0}-{\mathbf{z}}^{\ast}\|_{U}$ by Equation (33) and the monotonicity of $(\|{\mathbf{z}}^{j}-{\mathbf{z}}^{\ast}\|_{U})_{j\in{\mathbf{N}}}$ (Proposition 9). Thus, $\|{\mathbf{z}}_{\lambda}-{\mathbf{x}}\|_{U}\leq\|{\mathbf{z}}^{0}-{\mathbf{z}}^{\ast}\|_{U}+\|{\mathbf{z}}^{\ast}-{\mathbf{x}}\|_{U}.$ Therefore, for all $\lambda\in(0,1/\alpha]$, we have

Note that the upper key term identities (Proposition 14) and the fundamental inequality (Proposition 12) continue to hold when ${\mathbf{z}}^{k+1}$ is replaced by ${\mathbf{z}}_{\lambda}$. Thus, in each of the cases below, we will minimize the fundamental inequality over all $\lambda\in C$.

Part 2: First choose $\widetilde{\lambda}\in C$ small enough that $\rho+\mu-{\varepsilon}/{(\beta\widetilde{\lambda})}\leq 0$. Now recall that Proposition 14 proves the following inequality: $\kappa_{u}^{k}(\lambda)\leq\left(\rho-{\varepsilon}/({\beta\lambda})\right)\|{\mathbf{z}}^{k+1}-{\mathbf{z}}^{k}\|^{2}+2\gamma\lambda{\mathbf{g}}({\mathbf{x}}_{\mathbf{g}}^{k})-2\gamma\lambda{\mathbf{g}}({\mathbf{x}}_{\mathbf{f}}^{k}).$ Thus, the fundamental inequality, the cosine rule, and the identity $C=(0,1/\alpha]$ show

Part 3: We prove the result in the case that ${\mathbf{f}}$ is Lipschitz because the other case is symmetric. Proposition 14 proves the following identity: $\kappa_{u}^{k}(\lambda)=\left(1-{2}/{\lambda}\right)\|{\mathbf{z}}_{\lambda}-{\mathbf{z}}^{k}\|_{U}^{2}.$ Thus, the fundamental inequality, the cosine rule, and the identities ${\mathbf{x}}_{\mathbf{f}}^{k}-{\mathbf{x}}_{\mathbf{g}}^{k}=(1/\lambda)({\mathbf{z}}_{\lambda}-{\mathbf{z}}^{k})=T{\mathbf{z}}^{k}-{\mathbf{z}}^{k}$, ${\mathbf{x}}_{\mathbf{g}}^{k}-{\mathbf{x}}_{\mathbf{S}}^{k}=w({\mathbf{x}}_{\mathbf{g}}^{k}-{\mathbf{x}}_{\mathbf{f}}^{k})$, and $C=(0,2]$ show

Note that we can immediately strengthen the convergence result for PRS in Theorems 18 and 16. Indeed, we only need to assume that ${\mathbf{f}}$ or ${\mathbf{g}}$ is Lipschitz on the closed ball $\overline{B_{U}({\mathbf{x}}^{\ast};\|{\mathbf{z}}^{0}-{\mathbf{z}}^{\ast}\|_{U})}$ (where ${\mathbf{x}}^{\ast}=J_{\gamma U^{-1}(\partial{\mathbf{g}}+(1-w){\mathbf{S}})}({\mathbf{z}}^{\ast})$) of radius $\|{\mathbf{z}}^{0}-{\mathbf{z}}^{\ast}\|_{U}$ (under the metric $\|\cdot\|_{U}$) because for all $k\in{\mathbf{N}}$,

and by a similar derivation, $\|{\mathbf{x}}_{\mathbf{f}}^{k}-{\mathbf{x}}^{\ast}\|_{U}\leq\|{\mathbf{z}}^{0}-{\mathbf{z}}^{\ast}\|_{U}$. Thus, the sequences lie in the ball: $({\mathbf{x}}_{\mathbf{f}}^{j})_{j\in{\mathbf{N}}},({\mathbf{x}}_{\mathbf{g}}^{j})_{j\in{\mathbf{N}}}\subseteq\overline{B_{U}({\mathbf{x}}^{\ast},\|{\mathbf{z}}^{0}-{\mathbf{z}}^{\ast}\|_{U})}$. We also have $(\overline{{\mathbf{x}}}_{\mathbf{f}}^{j})_{j\in{\mathbf{N}}},(\overline{{\mathbf{x}}}_{\mathbf{g}}^{j})_{j\in{\mathbf{N}}}\subseteq\overline{B_{U}({\mathbf{x}}^{\ast},\|{\mathbf{z}}^{0}-{\mathbf{z}}^{\ast}\|_{U})}$ by the convexity of the ball. See [2, Proposition 8.28] for conditions that ensure Lipschitz continuity of convex functions on balls.

In general, the $o(1/\sqrt{k+1})$ convergence rates in Theorem 18 are the best PRS can obtain for $({{\mathbf{x}}}_{\mathbf{g}}^{j})_{j\in{\mathbf{N}}}$ [24, Theorem 11].

Applications

In this section we will show that the four algorithms from Section 2.2 are capable of solving highly structured optimization problems:

Let ${\mathcal{H}}_{0}$ be a Hilbert space, and let $f,g:\Gamma_{0}({\mathcal{H}}_{0})$. Let $n\in{\mathbf{N}}\backslash\{0\}$, and for $i=1,\cdots,n$, let ${\mathcal{H}}_{i}$ be a Hilbert space, let $h_{i},l_{i}\in\Gamma_{0}({\mathcal{H}}_{i})$, suppose that $h_{i}\square l_{i}\in\Gamma_{0}({\mathcal{H}}_{i})$, and let $B_{i}:{\mathcal{H}}_{0}\rightarrow{\mathcal{H}}_{i}$ be a bounded linear map. Finally, let ${\mathbf{B}}:{\mathcal{H}}_{0}\rightarrow\prod_{i=1}^{n}{\mathcal{H}}_{i}$ be the map $x\mapsto(B_{1}x,\cdots,B_{n}x)$. Then our model problem is as follows: