A first-order primal-dual method with adaptivity to local smoothness
Year: 2021

Abstract:
We consider the problem of finding a saddle point for the convex-concave objective minxmaxy f(x) + 〈Ax, y〉 − g∗(y), where f is a convex function with locally Lipschitz gradient and g is convex and possibly non-smooth. We propose an adaptive version of the Condat-Vũ algorithm, which alternates between primal gradient steps and dual proximal steps. The method achieves stepsize adaptivity through a simple rule involving ‖A‖ and the norm of recently computed gradients of f . Under standard assumptions, we prove an O(k−1) ergodic convergence rate. Furthermore, when f is also locally strongly convex and A has full row rank we show that our method converges with a linear rate. Numerical experiments are provided for illustrating the practical performance of the algorithm.

heading: 1 Introduction
text: In this paper we study a particular instance of the composite minimization problem
min x∈X f(x) + g(Ax), (1)
where f and g are convex, proper and lower-semicontinuous (l.s.c.), and A is a linear operator.
Problems of the form (1) have been studied in the literature under various assumptions on f and g. For the particular instances where g ◦ A is proximal-friendly and f is L-smooth, the objective is suitable for applying forward-backward splitting algorithms like the Proximal Gradient algorithm and its accelerated counterpart [Nesterov, 2013, Beck and Teboulle, 2009]. In general, however, the proximal operator of g ◦A is not easily computable and in such cases a popular approach is to decouple A and g by reformulating problem (1) as a convex-concave saddle-point problem:
min x∈X max y∈Y 〈Ax, y〉+ f(x)− g∗(y), (2)
where g∗ denotes the Fenchel conjugate of g. Objective (2) is typically addressed by primal-dual splitting algorithms which, under strong duality, can recover the solution to the original problem (1). In the particular case when f and g are proximal-friendly and possibly non-smooth, a very popular method is the Primal-Dual Hybrid Gradient proposed in [Chambolle and Pock, 2011], which was further extended to handle an additional L-smooth component with the Condat-Vũ algorithm [Condat, 2013, Vũ, 2013]. Convergence rates for the latter are studied in [Chambolle and Pock, 2016a].
Together, these classes of algorithms cover a broad range of problems in diverse fields such as signal processing, machine learning, inverse problems, telecommunications and many others. As a result, a great amount of research effort has gone into addressing practical concerns such as robustness
35th Conference on Neural Information Processing Systems (NeurIPS 2021).
to inexact oracles, acceleration and automation of stepsize selection. For a comprehensive list of examples and theoretical details we refer the reader to review papers [Combettes and Pesquet, 2011, Parikh and Boyd, 2014, Komodakis and Pesquet, 2015, Chambolle and Pock, 2016b]. In this work, we focus on the line of investigation studying stepsize regime automation for primal-dual algorithms targeting problem (2).
In their basic form, primal-dual methods require as input stepsize parameters belonging to a designated interval of stability, which depends on problem specific constants like the global smoothness parameter L and ‖A ‖. Dependence on such constants is undesirable because they may be costly to compute and oftentimes one can only access upper-bound estimates, thus leading to overly-conservative stepsizes and slower convergence. Moreover, the need to know L for setting the stepsizes prevents these methods from being applied to functions which are not globally smooth.
Consequently, recent efforts have gone towards devising methods with adaptive stepsizes [Goldstein et al., 2013, 2015, Malitsky and Pock, 2018, Pedregosa and Gidel, 2018]. These approaches resort to linesearch for finding good stepsizes at every iteration, and exhibit improved practical performance. It thus appears that better convergence comes at the cost of an indeterminate number of extra steps spent in subprocedures aimed at finding appropriate stepsizes.
In this work, we study problem (2) under the assumption that∇f is locally Lipschitz continuous and g is proximal-friendly. To illustrate the motivation of our framework, we take a prototypical example in image processing:
min x
1 2 ‖Kx− b ‖2 + λ ‖Dx ‖2,1 , K ∈ Rm×d, D ∈ R2d×d,
where x is an image, K is a problem-specific measurement operator, b is the vector of (possibly noisy) observations and D is the discrete gradient operator and the regularization term represents the isotropic TV norm. In order to apply any of the aforementioned primal-dual algorithms, one needs to first choose how to decouple the linear operators. There are three options: decoupling with respect to K leaves us with having to compute the proximal operator of the TV norm for the primal step, which is an iterative procedure [Chambolle, 2004]. Decoupling D implies performing gradient steps on f , since in general its proximal operator is not efficient. Finally, decoupling with respect to both implies increasing the dimensionality of the dual variable to m+ 2d, which is problematic for large d and m. The sensible choice is the second one (i.e., decoupling D), and the question we seek to answer with this work is:
Does there exist a method for solving (2) that adapts to the local problem geometry without resorting to linesearch?
Our contribution is to propose a first-order primal-dual scheme that answers this question in the affirmative and is accompanied by theoretical convergence guarantees. Using standard analysis techniques we show an ergodic convergence of O(k−1) when ∇f is locally Lipschitz and g is proximal-friendly, and a linear convergence rate for the case when f is in addition locally strongly convex and A has full row rank. We provide numerical experiments for sparse logistic regression and image inpainting, as well as use our method as a heuristic for TV-regularized nonconvex phase retrieval.
The rest of the paper is structured as follows: Section 2 provides details about related work; Section 3 introduces notation, along with technical preliminaries and assumptions to be used in our analysis; Section 4 reports the main theoretical results alongside partial proofs; finally, partial numerical results are provided in Section 5 with the rest being deferred to the appendix due to lack of space.

heading: 2 Related Work
text: Adaptive Gradient Descent (GD) methods. Arguably the most widespread of optimization methods, GD presents similar shortcomings for setting the stepsize as those described in the previous section. In particular, much research effort has gone in devising variants of the algorithm that remove the need to estimate the global smoothness constant L. In a recent work, Malitsky and Mishchenko [2020] propose an extremely simple and effective alternative for setting the stepsize τk adaptively at every iteration, as follows:
τk = min { τk−1 √ 1 +
τk−1 τk−2 , ‖xk − xk−1 ‖ 2 ‖∇f(xk)−∇f(xk−1) ‖
} . (3)
Adaptivity essentially comes ‘for free’ in (3), as it involves solely quantities which have already been computed. Moreover, the method requires only the weaker assumption of local smoothness, thus extending the reach of provably-convergent GD to a wider class of differentiable functions while maintaining the standard O(k−1) convergence rate. In this work we show that the above technique can be extended to the analysis of primal-dual methods, where it gives rise to an algorithm whose stepsizes adapt to the local geometry of the objective’s (locally) smooth component f .
Adaptive monotone variational inequality (VI) methods. Malitsky [2020] proposes an algorithm for solving monotone VIs with a stepsize that adapts to local smoothness similarly to (3). This method solves the very general formulation of finding u∗ such that 〈F (u∗), u−u∗〉+h(u)−h(u∗) ≥ 0, ∀u for a given monotone operator F which is locally Lipschitz continuous. Our template (2) can be recovered from theirs by setting u = (x, y), with
F (u) = F (x, y) = [ ∇f(x) +AT y −Ax ] ,
and h(u) = g∗(y). The advantages of this approach are the relaxed requirement of local Lipschitz continuity for F and the fact that knowledge of ‖A‖ is not required. However, since the VI framework is very general and does not take advantage of the problem structure (e.g. the fact that 〈Ax, y〉 is a bilinear term), the method comes with worse convergence bounds than algorithms specifically designed to solve (2). In addition, the algorithm requires as input an upper bound on the stepsizes, despite them being set in accordance to the estimated local smoothness.
First order primal dual algorithms and adaptive versions. A popular method for solving (2) when f is L-smooth is the Condat-Vũ algorithm (CVA) [Condat, 2013, Vũ, 2013]. The method’s convergence is subject to a global stepsize validity condition given by ( 1 τ − L ) 1 σ ≥ ‖A ‖
2, where τ and σ are the primal and dual stepsizes, respectively.
Another approach to solving problem (2) is via the Primal–Dual Fixed-Point algorithm based on the Proximity Operator (PDFP2O) or the Proximal Alternating Predictor–Corrector (PAPC) methods [Loris and Verhoeven, 2011, Chen et al., 2013, Drori et al., 2015]. This approach comes with less restrictive stepsize conditions than CVA owing to a different iteration style, but which nevertheless depend on the global smoothness constant L and ‖A ‖ and have to be carefully chosen. In order to alleviate the burden of choosing the stepsize parameters in CVA, Malitsky and Pock [2018] propose a linesearch procedure involving only dual variable updates and which, for certain problems such as regularized least squares, does not require any additional matrix-vector multiplications. A characteristic of this algorithm is that it maintains a constant ratio between primal and dual stepsizes through a hyperparameter β — a setup which we also use in this work.

heading: 3 Preliminaries
text: Consider problem (2) and let X ,Y be finite dimensional real vector spaces equipped with the standard inner product 〈·, ·〉 and the associated Euclidean norm ‖ · ‖ = √ 〈·, ·〉. We denote by g∗ the Fenchel conjugate of g in (1) defined as g∗(y) = supx{〈x, y〉− g(x)}. In order to not overload the ∗ notation, we use AT to denote the adjoint operator of A.
One can easily see that (2) is a primal-dual formulation of the following primal and dual optimization problems, of which the former is the same as (1):
min x∈X f(x) + g(Ax), max y∈Y −(f∗(−AT y) + g∗(y)).
A saddle-point (x∗, y∗) ∈ X × Y of problem (2) satisfies the following optimality conditions:
−AT y∗ = ∇f(x∗), Ax∗ ∈ ∂g∗(y∗). (4)
For (x′, y′) ∈ X × Y we define the following quantities:
Px′,y′(x) := f(x)− f(x′) + 〈x− x′, AT y′〉,
Dx′,y′(y) := g ∗(y)− g∗(y′)− 〈Ax′, y − y′〉,
Gx′,y′(x, y) := Px′,y′(x) +Dx′,y′(y). These functions are convex for fixed (x′, y′) and whenever (x′, y′) = (x∗, y∗), it holds that Px∗,y∗(x) ≥ 0, Dx∗,y∗(x) ≥ 0 and Gx∗,y∗(x, y) ≥ 0, with the latter quantity representing the primal-dual gap. We also define the gap restricted to a bounded subset B1 ×B2 ⊂ X × Y as:
GB1×B2(x, y) := sup (x′,y′)∈B1×B2 Px′,y′(x) +Dx′,y′(y),
and note that it is non-negative whenever B1 ×B2 contains a saddle-point. Given a function f : X → R and L > 0, we say that f is L-smooth if its gradient ∇f is Lipschitz continuous: ‖∇f(x) − ∇f(y)‖ ≤ L‖x− y‖,∀x, y. Furthermore, f is locally smooth if ∇f is Lipschitz continuous on any compact subset C: ∀C ⊂ X , ∃LC > 0 such that ‖∇f(x)−∇f(y)‖ ≤ LC‖x− y‖,∀x, y ∈ C. We also say that f is µ-strongly convex if f(y) ≥ f(x) + 〈∇f(x), y − x〉 + µ2 ‖x− y ‖ 2 ,∀x, y. Similarly, f is locally strongly convex if it is strongly convex on any compact subset C: ∀C ⊂ X , ∃µC > 0 such that f(y) ≥ f(x) + 〈∇f(x), y − x〉+ µC2 ‖x− y ‖ 2 ,∀x, y ∈ C. We define the proximal operator of a convex function g : X → R ∪ {∞} as proxg(x) = argminz { g(z) + 12 ‖x− z ‖ 2 }
, and say that g is ‘proximal-friendly’ if proxg(x) has a closed form solution or can be efficiently computed to high accuracy.
Finally, the following two blanket assumptions will hold throughout the paper: Assumption 3.1. Function f is convex and locally smooth, while g convex, l.s.c., and proximalfriendly. Assumption 3.2. A saddle-point exists for problem (2) and thus strong duality holds.
We note that Assumption 3.2 is standard in the literature (see e.g., [Chambolle and Pock, 2011]). Assumption 3.1, on the other hand, is weaker than the usual global L-smoothness premise and thus enlarges the category of admissible functions f with instances such as x 7→ exp(x). To illustrate, consider the aforementioned function defined on the reals: the global smoothness assumption clearly does not hold, however for any fixed interval [a, b] ⊂ R the smoothness constant can be chosen as exp(b).
For showing linear convergence of our method, we will add the following assumption: Assumption 3.3. Function f is locally strongly convex and operator A has full row-rank.

heading: 4 Algorithm and convergence
text: The primal-dual method proposed for solving problem (2) under assumptions 3.1 and 3.2, is provided in Algorithm 1 under the abbreviation APDA, which we use from here onwards. APDA follows the same structure as the basic CVA [Chambolle and Pock, 2016a] for the given assumptions. Notice that if we restrict Assumption 3.1 to L-smooth functions f , we can in fact recover CVA by setting θk = θ = 1 and τk = τ , σk = σ fixed such that ( 1 τ − L ) 1 σ ≥ ‖A ‖ 2.
Algorithm 1 Adaptive Primal Dual Algorithm (APDA) Input: x0 ∈ X , y0 ∈ Y, τinit > 0, τ0 =∞, θ0 = 1, β > 0, c ∈ (0, 1) x1 = x0 − τinit(∇f(x0) +AT y0) for k = 1, 2, . . . do
Set τk = min {
1 2 √ L2k+(β/(1−c))‖A ‖ 2 , τk−1
√ 1 + θk−1 } , σk = βτk, θk = τkτk−1
x̃k = xk + θk(xk − xk−1) yk+1 = proxσkg∗(yk + σkAx̃k)
xk+1 = xk − τk(∇f(xk) +AT yk+1) end for

heading: 4.1 High level ideas
text: We can rephrase the global stepsize condition of CVA by introducing a free parameter β > 0 which represents the ratio between the fixed dual and primal stepsizes: β = στ . With this change of variables, the stepsize validity condition becomes τ ∈ ( 0, 2
L+ √ L2+4β‖A ‖2
) .
Our algorithm disposes of CVA’s global condition and relies instead on a very similar but local criterion given by τk ∈ ( 0, 1
Lk+ √ L2k+2β‖A ‖ 2
) , where Lk :=
‖∇f(xk)−∇f(xk−1) ‖ ‖ xk−xk−1 ‖ provides an
estimate of the local smoothness constant and β = σkτk . In particular, this requirement is satisfied by the first part of the expression defining τk in APDA:
τk = min  1 2 √ L2k + (β/(1− c)) ‖A ‖ 2 , τk−1 √ 1 + θk−1  (5) where c ∈ (0, 1). Intuitively, this rule demands that τk does not overstep a constant related to the local curvature, thus allowing for larger stepsizes in flatter regions and correspondingly smaller ones otherwise.
By itself, the first term of (5) does not ensure convergence, since overly-aggressive and possibly destabilizing stepsizes might occur in near-linear regions. This issue is addressed by the second part of the expression (5) which, informally, prevents the stepsize from increasing ‘too fast’ in consecutive iterations. Specifically, the increase factor is at most √ 1 + θk−1, where θk =
τk−1 τk−2 .
Under these two local stepsize conditions we are able to show APDA’s convergence using the weaker assumption of local smoothness of f , thus conveniently removing the need of estimating a global smoothness constant L. Remark 4.1. While τk does not adapt to ‖A ‖, for many practical problems this fact is not a big hindrance. Function f typically represents the data fidelity term, whose smoothness constant L (should it exist) can far exceed ‖A ‖ – the linear operator enforcing structured regularization on x. A specific example are TV-regularized imaging problems, where A is the discrete gradient operator whose norm is bounded by √ 8 [Chambolle, 2004], while the data fidelity term may involve a very large number of measurements and a larger norm, consequently. Remark 4.2. APDA takes an additional primal step prior to the for-loop, which is controlled by τinit given as input. This is needed for estimating L1 in the first iteration. In practice we set τinit = 1e-9, a sufficiently small value to ensure that x1 does not depart too far from x0 and yield a good estimate of L1. Furthermore, the setting of τ0 = ∞ simply ensures that in the first step, τ1 = 1
2 √ L21+(β/(1−c))‖A ‖ 2 and has no impact on further steps. Finally, in our experiments we set
c = 1e-15 – this is a parameter introduced for theoretical purposes as explained in the following section.

heading: 4.2 Analysis – the base case
text: In short, the main steps of our analysis are: first, we establish the inequality that characterizes the dynamics of APDA given in Lemma 4.1 below. Based on it, we are able to prove the boundedness of sequences {xk} and {yk} in Theorem 4.1. In turn, sequence boundedness alongside the local smoothness property of f allows us to conclude that there exists a constant L > 0 such that f is L-smooth on the compact set Conv({x∗, x0, x1, . . .}) – the closed convex hull generated by {x∗, x0, x1, . . .}. Finally, we leverage this information to show that (xk, yk) converges to a saddle point of (2) and derive the associated ergodic convergence rates presented in Theorem 4.1. Lemma 4.1. Consider APDA along with Assumptions 3.1 and 3.2 and (x, y) ∈ X ×Y . Then, for all k and ηk ∈ ( βτk‖A ‖
1−c , 1−2τkLk 2τk‖A ‖
) ,
‖xk+1 − x ‖2 + 1
β ‖ yk+1 − y ‖2 + (1− ηkτk ‖A ‖ − τkLk) ‖xk+1 − xk ‖2
+ ηk − τkβ ‖A ‖
βηk ‖ yk+1 − yk ‖2 + 2τk(1 + θk)Px,y(xk) + 2τkDx,y(yk+1)
≤ ‖xk − x ‖2 + 1
β ‖ yk − y ‖2 + τkLk ‖xk − xk−1 ‖2 + 2τkθkPx,y(xk−1).
Moreover, it holds that:
1) τkLk < 12 < 1− ηkτk ‖A ‖ − τkLk,
2) 1β − τk‖A ‖ ηk > cβ > 0.
Proof sketch. The full proof is deferred to the appendix. We use algebraic manipulations, APDA’s update rules, the Cauchy-Schwarz and Young inequalities and properties of the prox operator to get the recurrence:
‖xk+1 − x ‖2 + 1
β ‖ yk+1 − y ‖2 + (1− τk ‖A ‖ ηk − τkLk) ‖xk+1 − xk ‖2
+
( 1
β − τk ‖A ‖
ηk
) ‖ yk+1 − yk ‖2 + 2τk(1 + θk)Px,y(xk) + 2τkDx,y(yk+1)
≤ ‖xk − x ‖2 + 1
β ‖ yk − y ‖2 + τkLk ‖xk − xk−1 ‖2 + 2τkθkPx,y(xk−1), (6)
where ηk > 0 is a free iteration-dependent constant involved in Young’s inequality.
In order to obtain anything worthwhile we would like to set ηk such that, when unrolling (6) over the iterations, the terms containing ‖xk+1 − xk ‖2 and ‖ yk+1 − yk ‖2 accumulate on the LHS with positive coefficients. More precisely, we ask that:{
1 β − τk‖A ‖ ηk > cβ , 1− τk ‖A ‖ ηk − τkLk > 12 , (7)
where c ∈ (0, 1). We note that the RHS of the first inequality could have been chosen as 0, however, we made it strictly positive due to technical reasons related to controlling the sequence ‖ yk+1 − yk ‖2. In practice, we choose c to be as small as possible.
A similar remark holds for the second inequality, where it would have been sufficient to set its RHS to τk+1Lk+1. Since this would considerably complicate the analysis, we make the observation that τkLk < 1 2 , ∀k and use this simpler uniform upper-bound instead.
The inequalities (7) are equivalent to asking that ηk ∈ ( τkβ‖A ‖
1−c , 1−2τkLk 2τk‖A ‖
) and what is left to show
is that this is a valid interval i.e., that the left endpoint is strictly smaller than its right counterpart. This condition amounts to solving a quadratic inequality in τk, whose solutions lie in the interval( 0, 1
Lk+ √ L2k+2(β/(1−c))‖A ‖ 2
) . The proof is concluded by showing that our choice of τk indeed
satisfies this constraint.
We are now ready to state the main convergence result in Theorem 4.1 below, whose full proof is given in the appendix.
Theorem 4.1. Consider APDA along with Assumptions 3.1 and 3.2, and let (x∗, y∗) ∈ X × Y be a saddle point of problem (2). Then, for all k
1) Boundedness. The sequence {(xk, yk)} is bounded. Specifically, for all k,
‖xk − x∗ ‖2 + ‖ yk − y∗ ‖2 ≤M,
where M := ‖x1 − x∗ ‖2 + 1β ‖ y1 − y∗ ‖ 2 + 12 ‖x1 − x0 ‖ 2 <∞.
2) Convergence to a saddle point. The sequence {(xk, yk)} converges to a saddle point of (2).
3) Ergodic convergence. Let Sk := k∑ i=1 τi, Xk := 1 Sk
( τk(1 + θk)xk +
k−1∑ i=1
(τi(1 + θi)− τi+1θi+1)xi ) and Yk := 1
Sk k∑ i=1 τiyi+1. Then, for any bounded
B1 ×B2 ∈ X × Y and for all k,
GB1×B2(Xk, Yk) ≤ M(B1, B2)
√ L2 + (β/(1− c)) ‖A ‖2
k ,
where L is the Lipschitz constant of ∇f over the compact set Conv({x∗, x0, x1, . . .}) and M(B1, B2) = sup(x,y)∈B1×B2 ‖x1 − x ‖ 2 + 1β ‖ y1 − x ‖ 2 + 12 ‖x1 − x0 ‖ 2.
The boundedness result of Theorem 4.1 point 1) implies that the closed set C = Conv({x∗, x0, x1, . . .}) is also bounded and hence compact. The local smoothness assumption on f then ensures that there exists L > 0 such that f is L-smooth over C. Note that such an L exists for any x0, y0 since the boundedness result itself holds for any initial conditions (though the value of such L cannot be generally known, as it is path-dependent). Using this fact, we can show a uniform
lower-bound on the primal stepsize: τk ≥ 1
2
( L2 + (β/(1− c)) ‖A ‖2 )−1/2 > 0, ∀k, which is
instrumental in deriving the subsequent convergence results, as well as Theorem 4.2. We emphasize that the appearance of constant L in the provided rates is a consequence of iterate boundedness, whose proof does not require its knowledge. Finally, we note that our rate is comparable to that of CVA in terms of constants.

heading: 4.3 Analysis under the additional Assumption 3.3
text: We now study APDA under the additional assumption of locally strongly convex f and full row rank A. Before proving the result of Theorem 4.2, a few remarks are in order. First, the boundedness result of Theorem 4.1 point 1) also holds for constant c = 0, since this constant was required only for proving convergence to a saddle point in point 2) of the theorem. Second, taking a smaller stepsize than the originally defined τk will not change the validity of Lemma 4.1 or the boundedness result of Theorem 4.1, as it remains within the required interval mentioned in section 4.1.
Consequently, for studying APDA under the additional Assumption 3.3 we can simplify the stepsize expression by taking c = 0, because now we are able to show iterate convergence directly by using the strong convexity and full row-rank assumptions. Specifically, we consider the stepsize:
τk = min  1 2 √ 4L2k + β ‖A ‖ 2 , τk−1 √ 1 + θk−1/2  , (8) which is smaller than the one originally considered and, due to the aforementioned remarks it ensures that APDA produces a bounded sequence. It follows that, under the local smoothness and local strong convexity assumptions, there exist constants L and µ such that f is L-smooth and µ-strongly convex over Conv({x∗, x0, x1, . . .}). The existence of these constants along with A being full row rank, in turn, allows us to derive a strengthened version of the inequality in Lemma 4.1 for (x, y) = (x∗, y∗):
‖xk+1 − x∗ ‖2 + ( 1
β + q1
) ‖ yk+1 − y∗ ‖2 + ( 1
2 + q2
) ‖xk − xk+1 ‖2 + q3 ‖ yk+1 − yk ‖2
+ 2τk(1 + θk)Px∗,y∗(xk) + 2τkDx∗,y∗(yk+1)
≤ (1− q4) ‖xk − x∗ ‖2 + 1
β ‖ yk − y∗ ‖2 +
( 1
2 − q5
) ‖xk − xk−1 ‖2 + 2τkθkPx∗,y∗(xk−1),
where q1, q2, q3, q4, q5 > 0 are constants given in the appendix. This new inequality represents in fact a contraction, which guarantees the linear convergence rate stated in Theorem 4.2.
Theorem 4.2. Consider APDA along with Assumptions 3.1, 3.2 and 3.3. Let (x∗, y∗) ∈ X × Y be a saddle point of problem (2). Furthermore, let τk be defined by (8) and let s := √ 4L2 + β ‖A ‖2
and t := √ 4µ2 + β ‖A ‖2, where µ, L are the strong convexity and smoothness constants of f over
the compact set Conv({x∗, x0, x1, . . .}).
Then, for all k:
‖xk − x∗ ‖2 + 1
β ‖ yk − y∗ ‖2 ≤ (1−min {p, q, r})kM,
where the rate constants are given by:
p = 1
2 , q =
µ 4s , r =
βσ2min(A)µ
βσ2min(A)µ+ 8s 2t+ 4L2s
,
and M = ‖x2 − x∗ ‖2 + ( 1 β + T ) ‖ y2 − y∗ ‖2 + 12 ‖x2 − x1 ‖ 2 + 2τ1Px∗,y∗(x1), T =
σ2min(A)µ
8s2t+ 4L2s , with σmin(A) representing the smallest singular value of A.
A few remarks are in order: first, as a sanity check, we observe that when A = 0 we recover the contraction factor of [Malitsky and Mishchenko, 2020] which is equal to q.
Second, we make some notes on how our rate compares with existing ones. To our knowledge, there are no explicit results regarding the linear convergence of CVA under assumptions similar to ours (linear rates are usually shown for the 3-component objective without assumptions on A — see e.g., [Chambolle and Pock, 2016a]). However, in the case of L-smooth and µ-strongly-convex f and full row-rank A, Chen et al. [2013] show the linear convergence of PDFP2O with rate:
‖xk − x∗‖2 ≤ ( ‖x1 − x0‖2 +
1
σmax(A) ‖y1 − y0‖2
)( 1−min { σ2min(A)
σ2max(A) , µ L
})k−1 ,
The rate presented in Theorem 4.2 has a comparatively worse contraction factor. The reason is that our iteration is set up in the style of CVA, where we essentially have a single stepsize to compute using the rephrasing from Section 4.1. Therefore, τk needs to obey the problem structure with respect to both L and ‖A‖, resulting in the ‘mixed’ term appearing in the denominator. Keeping the above in mind, the interested reader may find in the appendix that constants q and r come from a product between τk and other condition number-related quantities, which is tightly liked to the structure of the main inequality used in the paper. This makes the nice separation of condition numbers achieved in PDFP2O’s rate not possible in our case and, it seems, the analysis necessary to achieve this kind of adaptivity comes at the cost of worse constants (the same remark holds for [Malitsky and Mishchenko, 2020]).
PDFP2O, on the other hand, achieves a clean bound by having a different iteration style than CVA, as well as a fundamentally different kind of analysis where the iteration is expressed in fixed-point form to show convergence. In this context the stability conditions on the stepsizes are also relaxed — specifically, 0 < λ ≤ 1/σ2max(A) and 0 < γ < 2L in [Chen et al., 2013]. A drawback of this approach, however, is that the algorithm has no rate guarantees when f is only smooth and not strongly convex and only asymptotic convergence is shown. Also, PDFP2O requires 3 matrix-vector multiplications per iteration whereas we only require 2.

heading: 5 Experiments
text: We now present some numerical experiments conducted for APDA1. Additional problems and results are included the appendix. The experiments were implemented in Python 3.9 and executed on a MacBook Pro with 32 GB RAM and a 2,9 GHz 6-Core Intel Core i9 processor.
1See https://github.com/mvladarean/adaptive_pda.
The baseline we compare against in this section as well as the appendix is CVA, for which we use Algorithm 1 in [Chambolle and Pock, 2016a] (using g ≡ 0). In the particular case of sparse logistic regression we also compare against FISTA [Beck and Teboulle, 2009]. For obtaining x∗ we ran one of the algorithms for a large number of iterations.
5.1 Sparse binary logistic regression
We consider the problem of sparse binary Logistic Regression on 4 LIBSVM datasets [Chang and Lin, 2011] and show that adaptivity provides faster convergence in 3 of these cases. The objective we
consider is:
min x∈Rd F (x) := m∑ i=1
log(1 + exp(−bi〈qi, x〉))︸ ︷︷ ︸ f +λ ‖x ‖1︸ ︷︷ ︸ g , (9)
where (qi, bi) ∈ Rd × {−1, 1} and λ is the regularization parameter. APDA and CVA can be applied to this problem by setting A = I in formulation (2). Primal-dual algorithms are not the typical choice for solving (9), which is usually addressed by methods such as Proximal Gradient or FISTA [Beck and Teboulle, 2009]. However, we note that the computational costs of APDA and FISTA are comparable since the matrix-vector multiplication cost of the former is removed due to a A = I .
We choose λ = 0.005 ∥∥QT b∥∥∞, where QT = [qT1 , . . . qTm]T . For APDA we perform a parameter sweep over β ∈ [1e-3, 1e6] for each dataset and settle for: β = 2.68e3 for ijcnn; β = 5.18e4 for a9a; β = 3.16e1 for mushrooms; β = 3.73e-1 for covtype.
For CVA we sweep p ∈ [1e-3, 1e6] and set τ = 1‖A ‖/p+L and σ = 1p‖A ‖ — by construction, these stepsizes satisfy the validity condition and are as large as possible since the condition is satisfied with equality. We do an additional tuning procedure where we choose constants τ ∈ [1e-10, 1e2] and ξ ∈ [1e-5, 1e2] and set σ = τξ, which are subject to verifying the stepsize validity condition of CVA. Finally we select the best stepsizes across the two tuning phases to be (truncated to 3 decimals): τ = 9.869e-4, σ = 1.125e1 for ijcnn; τ = 2.655e-4, σ = 7.896e1 for a9a; τ = 9.936e-4, σ = 5.878e0 for mushrooms; τ = 7.728e-06, σ = 1e-06 for covtype.
Note that the Hessian of f is given by ∇2f(x) = QTD(x)Q, where D(x) is a diagonal matrix such thatDi,i(x) = σi(x)(1−σi(x)), where σi(x) = 11+exp(−bi〈qi,x〉) ∈ (0, 1). Clearly, over any compact set in C ⊂ X there exist Dmin := mini,x∈C Di,i(x) ∈ (0, 1) such that DminQTQ QTD(x)Q. As a result, a sufficient condition for local strong convexity is that the minimum eigenvalue of QTQ be greater than 0.
The convergence results along with stepsize comparison plots are presented in Figure 1. For dataset ijcnn we run APDA with the modified τk used in Theorem 4.2, since λmin(QTQ) = 75.13 and A has full rank. In the latter case, the legend identifier is APDA-strcnv. For the remaining datasets we use only the basic setting for τk, as λmin(QTQ) ≤ 1e-13. While APDA outperforms FISTA and CVA on ijcnn, a9a and mushrooms, it shows a relatively poor performance on covtype. We hypothesize that this is related to the condition number of QTQ, which is almost three orders of magnitude larger in the latter case: 9.2e22 versus 5.3e1, 2e20 and 2e17 for ijcnn, mushrooms and a9a, respectively. A similar behavior is seen in Figure 1.(c) of [Malitsky and Mishchenko, 2020].
Finally, the adaptive property of APDA’s stepsizes is visible in the stepsize comparison plots where they are shown to oscillate within at least one order of magnitude throughout the optimization process.

heading: 6 Limitations of APDA
text: The experiments presented in this paper (Section 5 and Appendix A) have the common trait of not imposing hard constraints on the primal variables. As a consequence, we are able to take plain gradient steps in the primal domain. However, for instances such as Poisson linear inverse problems [Bertero et al., 2009], the iterates xk need to reside in Rn+ because the primal objective contains log functions. APDA cannot handle such cases, as any constraints imposed on the primal variables will only be satisfied asymptotically. We consider such scenarios as a future research direction.

heading: Acknowledgements
text: The first author is grateful to Ya-Ping Hsieh for his feedback on the manuscript and for helpful research discussions throughout the development of this work. The authors also sincerely thank the anonymous reviewers for their time and their thoughtful, constructive feedback which helped improve and clarify this manuscript.
This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement no. 725594 - time-data), the Wallenberg Al, Autonomous Systems and Software Program (WASP) funded by the Knut and Alice Wallenberg Foundation, with the project number 305286. The work was also sponsored by the Department of the Navy, Office of Naval Research (ONR) under a grant number N62909-17-1-2111; by the Army Research Office and was accomplished under Grant Number W911NF-19-1-0404; by the Hasler Foundation Program: Cyber Human Systems (project number 16066). This work was also supported by the Swiss National Science Foundation (SNSF) under grant number 200021_178865/1.

Review summaries:
This paper studies a saddle-point problem with bilinear coupling and no constraints on primal as well as dual variables. Due to the unboundedness of the underlying problem, the paper assumes local Lipschitz smoothness for any compact sets. It proposes an adaptive version of the primal-dual method that estimates the local Lipschitz smoothness constant of the primal function. The authors show that the iterates remain in a bounded set around the saddle-point. Consequently, they can show asymptotic convergence of the limit point. They also show a convergence rate for a weaker notion of the primal-dual gap.

The authors propose an adaptive stepsize updating rule for the primal-dual Chambolle-Pock method that solves min f(x)+g(Ax) for locally-smooth f and cheap-prox-step g. The updating rule for adaptive stepsizes is simple in the sense that it only uses local gradients and ||A||, without needing any further subroutines like line search. They provide theoretical convergence rate for their proposed method under the local-smooth assumption and also with further local-strongly-convex assumption. They also corroborated their findings with experiment results that shows favorable  behavior of the new step-sizes.

This paper studies adaptivity of first-order primal-dual methods with respect to the smoothness or local smoothness of the problem. The manuscript builds upon a previous publication (Yura Malitsky and Konstantin Mishchenko. Adaptive gradient descent without descent. In Interna347 tional Conference on Machine Learning, pages 6702–6712. PMLR, 2020) to primal dual setup. Ergodic convergence is shown for the proposed method, where no knowledge of L is needed. The results is based on local smoothness which is defined for pairs of iterates of the sequence of points generated by the algorithm. The case of strongly convex function is also studied.