Accelerated Quasi-Newton Proximal Extragradient: Faster Rate for Smooth Convex Optimization
Authors: Ruichen Jiang
Aryan Mokhtari
Venue: NeurIPS 2023 spotlight (2023)
Keywords: convex optimization
quasi-Newton methods
Monteiro-Svaiter acceleration
Nesterov's accelerated gradient
online learning
TL;DR: We propose a quasi-Newton-type method that converges faster than Nesterov's accelerated gradient in the "small dimension" regime, while matching its optimal rate in the "large dimension" regime.

Abstract: In this paper, we propose an accelerated quasi-Newton proximal extragradient method for solving unconstrained smooth convex optimization problems. With access only to the gradients of the objective, we prove that our method can achieve a convergence rate of $\mathcal{O}\bigl(\min\\{\frac{1}{k^2}, \frac{\sqrt{d\log k}}{k^{2.5}}\\}\bigr)$, where $d$ is the problem dimension and $k$ is the number of iterations. In particular, in the regime where $k = \mathcal{O}(d)$, our method matches the _optimal rate_ of $\mathcal{O}(\frac{1}{k^2})$ by Nesterov's accelerated gradient (NAG). Moreover, in the the regime where $k = \Omega(d \log d)$, it outperforms NAG and converges at a _faster rate_ of $\mathcal{O}\bigl(\frac{\sqrt{d\log k}}{k^{2.5}}\bigr)$. To the best of our knowledge, this result is the first to demonstrate a provable gain for a quasi-Newton-type method over NAG in the convex setting.  To achieve such results, we build our method on a recent variant of the Monteiro-Svaiter acceleration framework and adopt an online learning perspective to update the Hessian approximation matrices, in which we relate the convergence rate of our method to the dynamic regret of a specific online convex optimization problem in the space of matrices.

Accelerated Quasi-Newton Proximal Extragradient: Faster Rate for Smooth Convex Optimization
Ruichen Jiang ECE Department The University of Texas Austin rjiangdutexas.er edu
Aryan Mokhtari ECE Department
The University of Texas at Austin mokhtari Cawttiinanni.n..m wtw.as.tee....mmmmmmmo.m edu
Abstract
In this paper, we present an accelerated quasi-Newton proximal extragradient method for solving unconstrained smooth convex optimization problems. With access only to the gradients of the objective function, we prove that our method
achieve convergence rate of O (min fz, vologk }). where d is the problem dimension and is the number of iterations. In particular, in the regime where = (d), our method matches the optimal rate of O f by Nesterov's accelerated gradient (NAG). Moreover, in the the regime where k = S2(d log d), it outperforms NAG and converges faster rate of O vdlosk) To the best of our knowledge, this result is the first to demonstrate a provable gain for a quasi-Newton-lyy method over NAG in the convex setting. To achieve such results, we build our method on a recent variant of the Monteiro-Svann adopt an online learning perspective to update the Hessian approximation matrices, in which we relate the convergence rate of our method to the dynamic regret of a specific online convex optimization problem in the space matrices.
acceleration framework and
Introduction
In this paper, we consider the following unconstrained convex minimization problem
min
f(x),
(1)
xERd
where the objective function f: Rd IR is convex and differentiable. We are particularly interested in quasi-Newton methods, which are among the most popular iterative methods for solving the problem in (I) (18 Like gradient descent and other first-order methods, quasi-Newton methods require only the objective's gradients update the iterates. On the other hand, they can better exploit the local curvature of by constructing Hessian approximation matrix and using it as a preconditioner, leading superior convergence performance. In particular, when the objective function in (I) is strictly convex strongly convex, it has long been proved that quasi-Newton methods achieve an asymptotic superlinear convergence rate (7115), which significantly improves the linear convergence rate obtained by first-order methods. More recently, there has been progress establishing local con-asymptotic superlinear rate of the form ((1/Vk)**) for classical quasi- Newton methods and their variants [16
However, all of the results above only apply under the restrictive assumption that the objective function f strictly or strongly convex. the more general setting where merely convex, the best of our knowledge, there is no result that demonstrates any form of convergence improvement by quasi-Newton methods over first-order methods. More precisely, accelerated gradient (NAG) 223 can achieve convergence rate of O(1/k²) if f is convex and has Lipschitz gradients. On the other hand, under the same setting, asymptotic convergence of
well known that Nesterov's
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

classical quasi-Newton methods has been shown in [3, 24] but no explicit rate is given. With certain conditions on the Hessian upproximation matrices, the works in 27 presented Sesiin-ettt--tppp methods with convergence rates of ( and O(1 /k²), respectively, which are no better than the rate of NAG and, in fact, can be even worse in terms of constants. This gap raises the following fundamental question: Can we design a quasi-Newton-Iype method that achieves a convergence rate faster than (1 for the smooth convex minimization problem in !
At first glance, this may seem impossible, as for any first-order method that has access only to a gradient oracle, one can construct a "worst-case" instance and establish lower bound of (011kk on the optimality gap 28 29]. It is worth noting that while such a lower bound is typically shown under a "linear span" assumption, i.e., the methods only query points in the span of the gradients they observe, this assumption is in fact not necessary and can be removed by the technique of resisting oracle (see, [30] Section 3.3]). In particular, this (2(1/²²) lower bound applies for any iterative method that only queries gradients of the objective, including quasi-Newton methods. On the other crucial assumption: it only works in the high-dimensional regime where the problem dimension d exceeds the number of iterations k. As such, it does not rule out the possibility of faster rate than (1//²²) when the number of iterations & is larger than d. Hence, ideally, we are looking for a method that attains the optimal rate of O11kk2) in the regime that k = (d) and surpasses this rate in the regime that l (dd)
hand, this lower bound is subject
Contributions. In this paper, we achieve the above goal by presenting an accelerated quasi-Newton proximal extragradient (A-QNPE) method. Specifically, under the assumptions that f in ① is convex and its gradient and Hessian are Lipschitz, we prove the following guarantees:
From any initialization, A-QNPE attain global convergence rate of O (minf 12, v 1). In particular, this implies that our method matches the optimal rate of O (FZZ) when (d), while converges at faster rate of O / when k- log d). Alternatively, we can bound the num- ber of iterations required to achieve an e-accurate solution by N, = O0ii( _ (log ب )2))).
In terms of computational cost, we show that the total number of gradient queries after N iterations can be bounded by 3.V. i.e., on average no more than per iteration. Moreover, the number of vector products to achieve an c-accurate solution can be bounded by (minn 20.s c0.225 }).
do
Combining the two results above, conclude that A-QNPE requires O Ominn ido. Li (log $ )2)) gradient queries to reach an e-accurate solution, which is at least as good as NAG and is further I )- To the best of our knowledge, this is the first result that demonstrates a provable advantage of -yatn--e method over NAG in terms of gradient oracle complexity in the smooth convex setting.
(log
superior when O 2 (og
To obtain these results, we significantly deviate from the classical quasi-Newtoe methods such as BFGS and DFP. Specifically, instead of mimicking Newton's 's method as in the classical updates, our A-QNPE method is built upon the celebrated Monteiro-Svaiter (MS) acceleration framework [3] 32)
which can be regarded as an inexact version of the accelerated proximal point method [3 34
Another major difference lies in the update rule of the Hessian approximation matrix. Classical quasi- Newton methods typically perform a low-rank update of the Hessian approximation matrix while enforcing the secant condition. On the contrary, our update rule is purely driven by our convergence analysis of the MS acceleration framework. In particular, inspired by (35, we assign certain loss functions to the Hessian approximation matrices and formulate the Hessian approximation matrix update as an online convex optimization problem in the space of matrices. Therefore, we propose to update the Hessian approximation matrices via an online learning algorithm.
Related work. The authors in (32) proposed a refined MS acceleration framework, which simplifies the line search subroutine in the original MS method [T.. By instantiating it with an adaptive second-order oracle, they presented an accelerated recond-order method that achieves the optimal rate of O (s.s). ). The framework in [32 serves as basis for our method, but we focus on the setting where we have access only to gradient oracle and we consider a Measi-Newtto-type update. Another (35], where the authors proposed a quasi-Newton proximal extragradient method with a global non-asymptotic superlinear rate. In particular, our Hessian approximation update is inspired by the online learning framework in 35) On the other hand, the major difference is that the authors in (35) focused on the case where f is strongly convex and presented global
closely related work
2

---

superlinear rate, while we consider the more general convex setting where f is only convex (may not be strongly convex). Moreover, we further incorporate the acceleration mechanism into our method, which greatly complicates the convergence analysis; see Remark 2 for more discussions.
Another class of optimization algorithms with better gradient oracle complexities than NAG cutting plane methods [36H44], which are distinct from quasi-Newton methods we study in this paper. In particular, in the regime where (Φ+) ), they can achieve the optimal gradient oracle complexity of Odloog ) [8 On the other hand, in the regime where = 2 /2 ), the complexity of the cutting plane methods is worse than NAG, while our proposed method matches the complexity of NAG.
2 Preliminaries
Next, we formally state the required assumptions for our main results.
Assumption 1. The function f is twice differentiable, convex, and Im-smooth. As result, we have 0 H fx) 3 for any x E Rd, where I E Rddd is the identity matrix.
Assumption 2. The Hessian of f is Lz-Lipschitz, i.e., we have || [| f(x) 2 f(y))oo S 121xx-ylb for any x, y E IRd, where |A _ ||px:||x||==1 ||Ax||2
We note that both assumptions are standard in the optimization literature and are satisfied by various loss functions such as the logistic loss and the log-sum-exy exp function (see, e.g., [6]].
Remark 1. We note that the additional assumption of Lipschitz Hessian does not alter the lower bound of (2(1/²2) that by a worst-case quadratic function, whose Hessian is constant (Lipschitz continuous with L2 = 0). Therefore, Assumption 2 does not eliminate this worst-case construction from the considered problem class, and thus the lower bound also applies to our setting.
discussed in the introduction. Indeed, this lower bound is established
Monteiro-Svann acceleration. As our proposed method uses ideas from the celebrated Monteiro- Svaiter (MS) acceleration algorithm [31], we first briefly recap this method. MS acceleration, also known as accelerated hybrid proximal extragradient (A-HPE), consists of intertwining sequences of iterates { ). } scalar variables } and {Aa} as well as step sizes [nk]. J. The algorithm has three main steps. In the first step, we compute the auxiliary iterate yı according to
Ok
7e
Th + 47
Y&
X +
Zk:
where
(2)
A
+ ak
A ak
In the second step, an inexact proximal point step x 2 Yk- Nk f fxx+11) is performed. To be precise, given parameter a é [0, 1), we find X++1 that satisfies
||x y& + T fx(x++1)) a o||x+11 - yull-
Then in the third step, the iterate z is updated by following the update
Zk+1 = Z ax fx(++)))
Finally, we update the scalar Ac+1 by Ax+1 A + ag- The above method has two implementation issues. First, to perform the update in (3) directly, one needs to solve the nonlinear system of equations x Yk + mr V fx() = 0 to a certain accuracy, which could be costly in general. To address this issue, a principled approach is to replace the gradient operator V fx) with a simpler approximation function ((x; yu) and select x as the (approximate) solution of the equation:
X - Yk + 1hs P(x =
use yk) V f(ys) and accordingly (4 is equivalent to x YR m Vf(y)), leading to the accelerated first-order method in BU. If we further have access to the Hessian oracle, we can use P(x; = ff + V²
For instance,
y&) and (4) becomes X leading to the second-ordes method in [3I.
Yk mh (I ny v2
However, approximating V fx) by P(x: leads to a second issue related to finding proper step size mh- More precisely, one needs to first select T/k, compute yı from (2), and then solve the system in 4 exactly approximately to obtain X However, these three variables, i.e., X ys and mn
may not satisfy the condition in S due to the gap between Vf(x) and P(x; that happens,
need to re- re-sct 77 and recalculate both yę and X until the condition in (3) is satisfied. To
3

---

address this issue, several bisection subroutines have been proposed in the literature 31 45 47 and they all incur a computational cost of lo/()/e) per iteration.
Optimal Monteiro-Svaiter acceleration. A recent paper [32 refines the MS acceleration algorithm by separating the update of Yı from the line search subroutine. In particular, in the first stage,
Te to compute ag and then Yı from (2. which will stay fixed throughout the line search scheme. In the second stage, we aim find a pair 8 and îix such that they satisty
| yı + i V >I of| yall-
To find that pair, we follow a similar line search scheme above, with the key difference that yı is fixed and ne can be different from 7/ that is used to compute yı. More precisely, for given fik. we find the solution of (4) denoted by X and check whether it satisfies (5) or not. If it does not, then we adapt the step size and redo the process until (5) is satisfied. Then given the values of these two parameters M7 and fu. the updates for and would change
describe next:
If nn 2 me, we update X++1 & AR+1 = A + oa and 2k+1 24 OK V (((+++)). Moreover, we increase the next tentative step size by choosing M]++1 = mk for some B E 10,
Otherwise, if fip < TJ.. the authors in 32 introduced a momentum damping mechanism. Define
Aa +744 & which is convex combination of x; and XK+1- Moreover, we update A&+1 = A + 7kaz and Zk+1 = Zk 7A. Vf Finally, we decrease the next tentative step size by choosing 7x+1 BT)-
(l- 7/k -ij//1 1. We then choose Xx+1 Actru Xx + TalAs
This approach not only simplifies the procedure by separating the update of {y&} } from the line search scheme, but it also shaves a factor of log(1/e) from the computational cost of the algorithm, leading to optimal first and second-order variants of the MS acceleration method. Therefore, as we will discuss in the next section, we build our method upon this more refined MS acceleration framework.
3 Accelerated Quasi-NNeton Proximal Extragradient
In this section, we present our accelerated quasi-Newton proximal extragradient (A-QNPE) method. An informal description of our method is provided in Algorithm I On a high level, our method can be Soaneirr-Svaiter-Neeee
method proposed in [32). In particular, we only query gradient oracle and choose the approximation function in (4) as P(x; yk) = f(yx) + B& yk). where B₃ is a Hessian approximation matrix obtained only using gradient information. Moreover, another central piece of our method is the update scheme of Bg. Instead of following the classical quasi-Newton updates such as BFGS or DFP, we use an online learning framework, where we choose sequence of matrices B₁ to achieve small dynamic regret for an online learning problem defined by our analysis; more details will be provided later in Section 3.2 We initialize our method by choosing Xo, Zo E Rd and setting Ao 0 and 170- do. where do is a user-specified parameter. Our method can be divided into the following four stages:
viewed as the quasi-Newton counterpart of the adaptive
In the first stage, we compute the scalar ag and the auxiliary iterate ya according to 2 using the step size Th- Note that Yı is then fixed throughout the k-th iteration.
In the second stage, given the Hessian approximation matrix B₃ and the iterate yk. we use a line search scheme to find the step size îjs and the iterate & such that
[&&+1 - yk + nI (V)(y)) + B₁ (x - yx))ll S on 011/ - ykll- | Ya + nM ( S (a) + 02)1 - yxll.
(6)
(7)
where a1 E 1).1) and a2 E (0, 1) are user-specified parameters with + a2 < 1. The first condition in (6) requires that &6+1 inexactly solves the linear system of equations (1+ Bx)
y&) + n V f(()) = O. where
E [0, 1) controls the accuracy. As a special case, we have xx+1 = yk II + the B VV f(ys) when a₁ O. The second condition in 7 directly comes
from (5) in the optimal MS acceleration framework, which ensures that we approximately follow the proximal point step & = yv iis V ( To find the pair (fn,, x satisfying both (6) and ( we implement a backtracking line search scheme. Specifically, for some 3 E (0, 1), we
iteratively try î/s = 3 for 2 and solve š from (6) until the condition in 7 satisfied.
The line search scheme will be discussed in more detail in Section|.11

---

Algorithm Accelerated Quasi-Newton Proximal Extragradient Method
1: Input: initial points Xo, 2o R initial step size au an, 02 E (0,1) with + < 1, (0.1)
2: Initialize: set Ao < 0 and Nu t do
3: for iteration & = O,.... N 1 do
4: Compute O2
Am
Ya + Ar Xx + Aat Zk
ak
5: Let be the largest possible step size in (ma'' : 0} such that
X1+1 Sai Ya - "VV( and || yv (f(*++1)) (a1+02)l yxll
6:
if fj 7/ then
7:
Set Rc+1, Z6+1 ua au (((e++))) Act1 Ak as
Set 7/411 6- fk/BB
9:
Set Bx+1 < Bs
10:
else
11:
Let 74 + na/TT
12:
Set X Actreak
(1-7elAg + 7ol (A Re+1, Ba+1 Ta (f((++1)....... Ae41 <- As + maa
13:
Set M <- 7/
14:
Set Wx f- [f((Ee) - (yr) ) and Sa t Re yA. where X is the last rejected iterate in LS
15:
Feed &a (B)
Beg
1
online learning algorithm and obtain B
17
16:
end if
17:
end for
In the third stage, we update the variables X Zk+1. A&+1 and set the step size 7 in the next iteration. Specifically, the update rule follow depends on the outcome of the line search scheme. Tik i.c., the line search scheme accepts the initial trial step size, we let
In the first case where ihe-
X &s+1s
Zk+1 Zμ a2 (8*++)) AR+1 = A + ak
(8)
as in the original MS acceleration framework. Moreover, this also suggests our choice of the step size Ti may be too conservative. Therefore, we increase the step size in the next iteration by 7h+11 f (.. In the second case where nn < 1J., i.e. the line search scheme backtracks, we adopt the momentum damping mechanism in 32
- Ak))A)A x
(e<AA + aa)
xx+1, Zk+1 = Y V f(((++))) Act1 = +
X
A4 + 7k.s
A + nak
(9)
where Tk = /m) < Accordingly, we decrease the step size in the next iteration by letting 1 1ix (note that ina < vR)).
In the fourth stage, we update the Hessian approximation matrix B&+1- Inspired by [35], we depart from the classical quasi-Newton methods and instead let the convergence analysis guide our update scheme. As will show in Section 3.2. the convergence rate of our method is closely related to the cumulative loss E (B,) incurred by our choices of (Bx), where B = {k: : î < 7 denotes the indices where the line search scheme backtracks. Moreover, the loss function has the form & (B²) |Wa_BB where w A fics) f((x). 5y - R2 - Xx and X₂ is an sull= auxiliary iterate returned by our line search scheme. Thus, this motivates us to employ an online learning algorithm to minimize the cumulative loss. Specifically, in the first case where î nkk (i.e., k & B3,, the current Hessian approximation matrix B, does not contribute to the cumulative loss and thus we keep unchanged (cf. Line 99. Otherwise, we follow an online learning algorithm in the space of matrices. The details will be discussed in Section 3.2
Finally, provide convergence result in the following Proposition for Algorithm I which serves as the basis for our convergence analysis. We note that the following results do not require additional conditions on B other than the ones in (6) and (7). The proof is available in Appendix A.
Proposition 1. Let {xx IN be the iterates generated by Algorithm If f is convex, we have
||zo -x112 and AN 2 -3-3
N
f(xN) - f(x')
2AN
4(2 - 2B)²
k=
5

---

Proposition| I characterizes the convergence rate of Algorithm[Iby the quantity AN. which can be further lower bounded in terms of the step sizes {. Moreover, we can observe that larger step sizes will lead to faster convergence rate. On the other hand, the step size rik is constrained by the condition in 77) which, in turn, depends on our choice of the Hessian approximation matrix B Thus, the central goal of our line search scheme and the Hessian approximation update is to make the step size fie as large as possible, which we will describe next.
3.1 Line Search Subroutine
In this section, we specify our line search subroutine to select the step size ih and the iterate Xx+1 in the second stage of A- A-QNPE. For simplicity, denote (yk) by g and drop the subscript in yı and B,. In light of (6) and 77. our goal in the second stage is to find a pair (f), x+) such that
y + ig + B(x y))ll _ a || yll,
(10)
| - y + / [*+)[ - (0r + o2)ll yll.
(11)
As mentioned in the previous section, the condition in (TO) can be satisfied by solving the linear system (I + (x+ y) = -ng to desired accuracy. Specifically, we let
8+ = LinearSolver(T + rjB - ig: mn) and =
(12)
where the oracle LinearSolver is defined as follows.
Definition 1. The oracle LinearSolver( A, b;) takes matrix A E $4 a vector b E Rd and a E (0, 1) as input, and returns an approximate solution S. satisfyyng ||As+ - bll S a||s+|]
The most direct way to implement LinearSolver(A. A, b; a) is to compute $ A A which however costs (dd)) arithmetic operations. Alternatively, we can implement the oracle more efficiently by using the conjugate residual method 48), which only requires computing matrix-vector products and thus incurs a cost of (dd). The details are discussed in Appendix E. We characterize the total number of required matrix vector products for this oracle in Theorem
Now we are ready to describe our line search scheme with the LinearSolver oracle (see also Subrou- tine in Appendix B). Specifically, we start with the step size 1i and then reduce it by a factor 3 until we find pair (f). x.) that satisfies (TT). can be shown that the line search scheme will terminate
in finite number of steps and return pair (f). x.) satisfying both conditions in (10) and (T1) (see Appendix B.I). Regarding the output, we distinguish two cases: (i) If we pass the test in (II) on our first attempt. we accept the initial step size and the corresponding iterate x+ (cf. Line 10 in Subroutine Db. (ii) Otherwise, along with the pair (f), &.). we also return an auxiliary iterate x, that we compute from 12) using the rejected step size n (cf. Line 12 in Subroutine [ As we shall see in Lemma] T the iterate X is used to derive a lower bound on 7/. which will be the key to our convergence analysis and guide our update of the Hessian approximation matrix. For ease of notation, let B be the set of iteration indices where the line search scheme backtracks, i.c., B {ke nn. < mp)-
Lemma 1. For k & B we have în TJk while for & B we have
02 yull
(1 + a1) 3(1 0n
f
and
| - 1oull 5
| yull-
( - Vf(x)) (x yx)ll
(13)
Lemma provides a lower bound on the step size ńn in terms of the approximation error (x*+++))
Vf(yx)
yl Hence, better Hessian approximation matrix B, leads to a larger step size, which in turn implics faster convergence. Also note that the lower bound uses the auxiliary iterate X&+1 that is not accepted as the actual iterate. Thus, the second inequality in (13) will be used to relate || - yull with | yx Il. Finally, we remark that to fully characterize the computational cost. we need to upper bound the total number of line search steps, each of which requires a call to LinearSolver and a call to the gradient oracle. This will be discussed in Theorem
Bx T[(yx)_BB(Xx..
3.2 Hessian Approximation Update via Online Learning with Dynamic Regret
In this section, we discuss how to update the Hessian approximation matrix B; in the fourth stage of A-QNPE. As mentioned earlier, instead of following the classical quasi-Newton updates, we directly motivate our update policy for Bü from the convergence analysis. The first step is to connect the

---

convergence rate of A A-QNPE with the Hessian approximation matrices (Ba}. By Proposition | we define the absolute constant C
if
_
2(2- then we can write
(1-v3)²
120 x=1
Cilla x1 N-1
Cilizo x21 N2.5
N-1
f(xN) - f(x){
(14)
2A
k=
R=
where the last inequality follows from Hölder's inequality. Furthermore, we can establish upper
N-1 1 in terms of the Hessian approximation matrices {Bx], show next. bound on kk=u oi
Lemma 2. Let {fn N-1
be the step sizes in Algorithm
using Subroutine
1 Then we have
k=0
N-1
2 32
32
I BeS 11
+ (1 - -- -32)0233
(15)
- 80)00
|8x 1
k=
U. USKS 1,KEB
where Wi
V f
V f(yx) and Sk.
x for E B.
The proof of Lemma 2 is given in Appendix C.
On high level, for those step sizes în with i E B, we can apply Lemma and directly obtain a lower bound in terms B On the other hand, for connect B). As result, note that the sum in (15) only involves the Hessian approximation matrices (B, : k B).
k & B, we have fi T/k and
update rule in Lines 8 and 13 of Algorithm aallows
the sequence {Nı Jk=
N-1 with the backtracked step sizes {i
In light of (14) and (15). our update for B₁ aims to make the right-hand side of (15) as small as possible. To achieve this, we adopt the online learning approach in (35) and view the sum ([5) as the cumulative loss incurred by our choice of [Bx]. To formalize, define the loss at iteration by
if & B, otherwise,
E B
W
Bag 12 1
(16)
and consider the following online learning problem: (i) At the k-th iteration, we choose B Z B L, (ii) We receive the loss function (B) defined in ( (iii) We update our Hessian approximation matrix to B&+1. Therefore, we propose to employ an online learning algorithm to update the Hessian approximation matrices {B&), and the task of proving a
where Z A (B E 54
convergence rate for our A A-QNPE algorithm boils down to analyzing the performance of our online
learning algorithm. In particular, upper bound on the cumulative loss En= (Bx) will directly translate into a convergence rate for A- A-QNPE by using [14 and (5).
Naturally, the first idea is to update B₁ by following projected online gradient descent 44 While this approach would indeed serve our purpose, its implementation could be computationally expensive. Specifically, like other projection-based methods, it requires computing the Euclidean projection onto the set Z each iteration, which our case amounts to performing a full d x d matrix eigendecomposition and would incur cost of (de) (see Appendix C.2 Inspired by the recent work in [50], we circumvent this issue by using a projection- free online learning algorithm, which relies on approximate separation oracle for Z instead of projection oracle. For simplicity, first translate and rescale the set Z via the transform B- 2 (B 41 I) to obtain 2 & (B Sd |B||op 1).
T.j
The approximate separation oracle SEP (P: 6, q) then defined as follows.
Definition 2.
The oracle SEP(W; 5. 9) takes a symmetric matrix W 65,8 8> and q (0, 1) as input and returns scalar > 0 and matrix S E Sd with one of the following possible outcomes:
Case I: S 1, which implies that, with probability at least 1 W E 2;
Case II: y > 1, which implies that, with probability at least 1 - W/r E 2, |S||e S 3 and
(S, w B) 2 7 - 1 for any B such that B 2.
To sum up. SEP(W; 6, q) has two possible outcomes: with probability 1 - either it certifies that W E, or it produces scaled version of W that belongs to å and an approximate separation hyperplane between W and the set 2. As show in Appendix E.2 implementing this oracle
7

---

requires computing the two extreme eigenvectors and eigenvalues of the matrix W inexactly, which can be implemented efficiently by the randomized Lanczos method [5T]
Building on the SEP(W; 6. oracle, we design a projection-free online learning algorithm adapted from [5 Subroutine 2]. Since the algorithm is similar the one proposed in [35), we relegate the
details to Appendix but sketch the main steps in the analysis. To upper bound the cumulative loss
N -N-1, k= & (Bx), we compare the performance of our online learning algorithm against sequence
of reference matrices {H&)
Specifically, aim control the dynamic regrer 49 52 53
defined by D-Reg ([H+) I Ex=O (4x()x) - Eg (H+)), as well as the the cumulative loss En N-1
&x(Hs) by the reference sequence. In particular, in our analysis we show that the choice of HA f(y()) for .,..., .V- 1 allows us to upper bound both quantities.
Remark While our online learning algorithm is similar to the one [35], our analysis more challenging due to the lack of strong convexity. Specifically, since I is assumed to be strongly convex (35. the iterates converge to x at least linearly, resulting in less variation the loss functions {
Hence, the authors in [35] let H₁ = H' fx*( for all & and proved that En l (H²) remains
bounded. In contrast, without linear convergence, we need to use a time-varying sequence } to
control the cumulative loss. This in turn requires us to bound the variation Σi=o |H&+1 - H IF.
which involves careful analysis of the stability property of the sequence {yk} in Algorithm
4 Complexity Analysis of A-QNPE
In this section, we present our main theoretical results: we establish the convergence rate of A-QNPE (Theorem D and characterize its computational cost in terms of gradient queries and matrix- product evaluations (Theorem The proofs are provided in Appendices D Dnd E.3
Theorem 1. Let } be the iterates generated by Algorithm using the line search scheme in Section 3.1 where 01,02 10,1) with on + 02 < and B (0, 1), and using the Hessian approximation update in Section 3.2 (the hyperparameters are given in Appendix D). Then with probability at least p. the following statements hold, where C₁ (i 4,...; 10) are absolute constants only depending on an, 02 and B.
(a) For any k 2 0, we have
f(xx) - f(x')= CaTn1120 1²
k2 + Czllzo - x²1
(17)
kg.5.
(b) Furthermore, for any k 0,
f(xx) - f(x)) 120 - 112
[max Ly
J1:2.5
M + LodLi - -x" logg
m3 op
k2..5
VM
(18)
where we define log* (ar) A max(log(x), 03 and the quantity M is given by
Co + +CL C BB f()(2)) + CoL x x² + CioLy - -x"1| x*1|. (19)
M
Both results in Theorem I are global, as they are valid for any initial points xo, Zo and any initial matrix Bo. Specifically, Part (a) of Theorem | shows that A- A-QNPE converges at rate of O(1 k22).
matching the rate of NAG 23 that is known to be optimal in the regime where k = O(d) 28 29
Furthermore, Part (b) of Theorem presents a convergence rate of O /oo((k))///// /12 2.5. To see this,
note that since we have 0- < By < I and 0 - f(z)) < I,I, in the worst case [ (z.))] in the expression of (9) can be upper bounded by d,d. Thus, assuming that L, Le and ||zg x"|
are on the order of O(1), we have M = (d) and the convergence rate in (18) can be simplified
to O( (( log (2.5). Notably, this rate surpasses the O(1/k²) rate when k = 2(d d). To the best of our knowledge, this is the first work to show convergence rate faster than (11k?²) for quasi-Newton-type method the convex setting, thus establishing provable advantage over NAG.
8

---

Table 1: The comparison of NAG and our proposed method in terms of computational cost.
Methous
Gradient queires
Matrix vector products
O( 0.9
A-QNPE (ours)
((min(e 0.6, 20:22 0.4))
(mmindd0.20 0.025))
Remark 3 (Iteration complexity). Based on Theorem [ we can find A-QNPE's iteration complexity.
Define N₂ as the number of iterations required by A A-QNPE to find an e-accurate solution, i.e., fxx - f(x') S When ve log the rate in (17 is better and have N₂ O ((,t5)
Conversely, when v allog the rate in (18) is the better one, resulting in N, O( log 2 )02).
Hence, to achieve an e-accurate solution A A-QNPE requires (minn 20.3.
di 20.4 (log _ ))))) iterations.
Remark (Special case). If the initial point Zo is close to an optimal solution x* and the initial Hessian approximation matrix Bg is chosen properly, the dependence on d in the convergence rate of (18) can be eliminated. Specifically, if l|zo x"1 = -(2) and we set B, = fzo). then we have M = (1) and this leads to a local dimension-
dimension-independden independent rate of O O(vlogk/k29)
Recall that in each iteration of Algorithm|1 [ we need to execute a line search subroutine (Section 3.1 and Hessian approximation update subroutine (Section 3.2p. Thus, to fully characterize the computational cost of Algorithm | we need to upper bound the total number of gradient queries as well as the total number of vector product evaluations, which is the goal of Theorem 2
Theorem 2. Recall that Ne denotes the minimum number of iterations required by Algorithum find an c-accurate solution according to Theorem 1 Then, with probability at least p.
1o
(a) The total number of gradient queries is bounded by 3N, + log1/a
CoLi
02
(b) The total number of matrix-vector product evaluations in the LinearSolver oracle is bounded
by Ne + CCnvaola + C1av L1
120-x'1 x''x'1 where C11 and C12 are absolute constants.
2
(c) The total number of matrix-vector product evaluations the SEP oracle is bounded by O (N! 25 (log N.) (og vāN.))
If the initial step size is chosen as ou 02 Theorem implies that A-QNPE requires more Z, than 3 gradient queries per iteration average. Thus. the gradient oracle complexity of A-QNPE
d the same as the iteration complexity, i.e., O(min a 4 (log d ).2))) On the other hand, the complexity in terms of meatrx-vvector products is worse.
More precisely, by using the expression of N Remark Parts (b) and (c) imply that the total number of matrix- vector product evaluations
LinearSolver and SEP oracles can be bounded by ( ((t.)) ) and Óminn 20.5 20.025 }), respectively. For casier comparison, summarize the detailed computational costs of NAG and our method A-QNPE to achieve an e-accuracy in Table T We observe that A-QNPE outperforms NAG in terms of gradient query complexity: It makes equal or fewer gradient queries especially when & da On the other hand, A A-QNPE requires additional matrix-vecte product computations to implement the LinearSolver and SEP oracles. While this is limitation of our method, in some cases, gradient evaluations are the main bottleneck and can be more expensive than matrix-vector products, As a concrete example, consider the finite-sum minimization problem f(x) 1 Σ fi(x). In this case, one gradient query typically costs (nd), while one matrix-vector product costs O(d)) Thus, the total computational cost of NAG and A-QNPE can be bounded by O ((15)) and Ó nd d2:25 ).
d0.25
20.4
respectively. In particular, our method incurs a lower computational cost when < 7 and n > dl.25_
5 Experiments
In this section, we compare the numerical performance of our proposed A A-QNPE method with NAG and the classical BFGS quasi-Newton method. For fair comparison, also use a line search

---

10
S0d TTT +- Bris wy
0.6
2 qUF
100
0.5
F36 F1x))
70
0000000
2000
0000
0.1
500
1050
1500
2000
2000
Gradlest
Gendbent queedes
(a) Convergence by iteration (b) Convergence by gradient queries
(c) Histogram of gradient queries
Figure 1: Numerical results for logistic regression on synthetic dataset.
10
Rolo SNPP
10
109
102
0.5
0.4
[(a)-
10
0.3
0.2
10
0.1
15
1.5
6
=100
quantes
x100
Gradicnt quecies
(a) Convergence by iteration
(h) Convergence by gradient queries
(c) Histogram of gradient queries
Figure 2: Numerical results for log-sum-exp function on a synthetic dataset.
We would like to highlight that our paper mainly focuses establishing provable gain for quasi-Newton methods with respect to NAG, and our experimental results are presented to numerically verify our theoretical findings. In the first experiment, we focus logistic regression problem with the loss function f(x) = E log(1 + e (a:,)), where a,...., a, E Rd are feature vectors and yI..... E {-1, are binary labels. We perform our numerical experiments on a synthetic dataset and the data generation process is described in Appendix F In the second experiment, we consider the
scheme NAG and BFGS obtain their best performance 54 15
log-sum-exp function fxx) = log(>> i= e( b.) ). where we generate the dataset {(a, b, }
(more details Appendix [F). As we observe in F [a) and Fig. Da), our proposed A A-QNPE method converges in much fewer iterations than NAG, while the best performance is achieved by BFGS. Due the use line search, we also compare these algorithms in terms of the total number of gradient queries. Moreover, additional plots in terms of the running time are included in Appendix F As illustrated in Fig. bb) and Fig. bb), A A-QNPE still outperforms NAG but the relative gain becomes less substantial. This is because the line search scheme in NAG only queries the function value at the new point, and thus it only requires one gradient per iteration. On the other hand, should add that the number of gradient queries per iteration for A-QNPE is still small as guaranteed by our theory. In particular, the histogram of gradient queries in Fig. Tec and Fig.2) ([c) shows that most of the iterations of A A-QNPE require 2-3 gradient queries with an average of less than 3. Finally, although there is no theoretical guarantee showing a convergence gain for BFGS with respect NAG, observe that BFGS outperforms all the other considered methods in experiments. Hence, studying the convergence behavior of BFGS (with line search) in the convex setting is an interesting research direction to explore.
following a similar procedure as in 16
6 Conclusions
We proposed a quasi-Newton variant of the accelerated proximal extragradient method for solving smooth convex optimization problems. We established two global convergence rates for our A-QNPE method, showing that it requires O(min 20s 2 }) gradient queries to find e-accurate solution. In particular, the regime where 2 dz ), A A-QPE achieves gradient oracle complexity of
d0.2
O (1))), ), matching the complexity of NAG. Morcover, in the regime where Ō ). it outperforms
NAG and improves the complexity to Ō((²)) ). To the best of our knowledge, this the first result showing provable gain for quasi-Newton-typee method over NAG in the convex setting.
10

---

Acknowledgments and Disclosure of Funding
The research of R. Jiang and A. Mokhtari is supported in part by NSF Grants 2007668, 2019844, and 2112471, ARO Grant W911NF2 10226, the Machine Learning Lab (MLL) at UT Austin, the Wireless Networking and Communications Group (WNCG) Industrial Affiliates Program, and the NSF AI Institute for Foundations of Machine Learning (IFML). The authors would also like to thank Yair Carmon and the anonymous reviewers for their comments on the first draft of the paper.
References
[1] W. C. Davidon. Variable metric method for minimization. Techinical Report ANL-5990, Argonne National Laboratory, Argonne, IL, 1959.
[2] Roger Fletcher and Michael JD Powell. A rapidly convergent descent method for minimization. The Computer Journal, 6(2):163-168 1963.
[3] Charles G Broyden. The convergence of single-rank quasi-newton methods. Mathematics of Computation, 24(1 14(133)):65-32,, 1970.
[4] Roger Fletcher. A new approach to variable metric algorithms. The computer journal, 13(3):317- 322, 1970.
[5] Donald Goldfarb. A family of variable-metric methods derived by variational means. Mathe- matics of computation, 24(109):23-26, 1970.
[6] David F Shanno. Conditioning of quasi-newton methods for function minimization. Mathemat- ics of computation, 24(111):6447-5
1970.
[7] Andrew R Conn, Nicholas IM Gould, and Ph Toint. Convergence of quasi-newton matrices generated by the symmetric rank one update. Mathematical programming, 50(1)1177-195, 1991.
[8] Fayez Khalfan, Richard H Byrd, and Robert B Schnabel. A theoretical and experimental study of the symmetric rank-one update. SIAM Journal on Optimization, 3(1):1-24, 1993.
[9] Charles George Broyden, John E Dennis Jr, and Jorge J Moré. On the local and superlinear convergence of quasi-newton methods. IMA Journal of Applied Mathematics, 12(3):22-2245 1973.
[10] John E Dennis and Jorge Moré. A characterization of superlinear convergence and its application to quasi-newtoe methods. Mathematics of computation, 28(126):549-560 1974.
[11] MJD Powell. On the convergence of the variable metric algorithm. IMA Journal of Applied Mathematics, 7(1):21-36, 1971.
[12] Laurence Charles Ward Dixon. Variable metric algorithms: necessary and sufficient condi- tions for identical behavior of nonquadratic functions. Journal of Optimization Theory and Applications, 10(1):34 40, 1972.
[13) Richard H Byrd, Jorge Nocedal, and Ya-Xiang Yuan. Global convergence of a class of quasi- newton methods on convex problems. SIAM Journal on Numerical Analysis, 2415):1 1987.
[14) Richard H Byrd, Humaid Fayez Khalfan, and Robert B Schnabel. Analysis of a symmetric rank-one trust region method. SIAM Journal on Optimization, 6(4):1025-1039, 1996.
[15] Jorge Nocedal and Stephen J. Wright. Numerical Optimization. Springer Science+Business Media, LLC, 2006.
[16] Anton Rodomanov and Yurii Nesterov. Greedy quasi-newton methods with explicit superlinear convergence. SIAM Journal on Optimization, (11)77885-11, 2021.
[17] Anton Rodomanov and Yurii Nesterov. Rates of superlinear convergence for classical quasi- newton methods. Mathematical Programming, 2021.

---

[18] Anton Rodomanov and Yurii Nesterov. New results on superlinear convergence of classical quasi-newton methods. Journal of Optimization Theory and Applications, 188())744-769, 2021.
[19] Qiujiang Jin and Aryan Mokhtari. Non-asymptotic superlinear convergence of standard quasi- newton methods. Mathematical Programming. 2022.
[20] Qiujiang Jin, Alec Koppel, Ketan Rajawat, and Aryan Mokhtari. Sharpened quasi-Newton methods: Faster superlinear rate and larger local convergence neighborhood. In Proceedings of the 39th International Conference on Machine Learning, volume 162 of Proceedings of Machine Learning Research, pages 10228-10250. PMLR, 2022.
[21] Dachao Lin, Haishan Ye, and Zhihua Zhang. Greedy and random quasi-newton methods with faster explicit superlinear convergence. Advances in Neural Information Processing Systems, 34:6646-6657 2021.
[22] Haishan Ye, Dachao Lin, Xiangyu Chang, and Zhihua Zhang. Towards explicit superlinear convergence rate for srl. Mathematical Programming, pages 1-31, 2022.
[23] Yurii Nesterov. A method of solving convex programming problem with convergence rate O(1/k2). Soviet Mathematics Dokludy, pages 372-376, 1983.
[24] M.J.D, Powell. Some properties of the variable metric algorithm. Numerical methods for nonlinear optimization, pages 1-17, 1972.
[25] Katya Scheinberg and Xiaocheng Tang. Practical inexact proximal quasi-newton method with global complexity analysis. Mathematical Programming, 160:495-529, 2016.
[26] Hiva Ghanbari and Katya Scheinberg. Proximal quasi-newton methods for regularized con- vex optimization with linear and accelerated sublinear convergence rates. Computational Optimization and Applications, 69:597-6211 2018.
[27] Dmitry Kamzolov, Klea Ziu, Artem Agafonov, and Martin Takáč, Accelerated adaptive cubic regularized quasi-newton methods. arXiv preprint (rXi::2302.04987,, 2023.
[28] A. Nemirovski and D Yudin. Problem Complexity and Method Efficiency in Optimization. Wiley, 1983.
[29] Yurii Nesterov. Lectures on Convex Optimization. Springer International Publishing, 2018.
[30] Yair Carmon, John C Duchi, Oliver Hinder, and Aaron Sidford. Lower bounds for finding stationary points i. Mathematical Programming. 184(1-2):11-100 2020.
[31] Renato DC Monteiro and Benar Fux Svaiter. An accelerated hybrid proximal extragradient method for convex optimization and its implications to second-order methods. SIAM Journal on Optimization, 1392) 2013.
[32] Yair Carmon, Danielle Hausler, Arun Jambulapati, Yujia Jin, and Aaron Sidford. Optimal and adaptive monteiro-svaiter acceleration. In Advances Neural Information Processing Systems, 2022.
[33] Osman Güler. New proximal point algorithms for convex minimization. SIAM Journal on Optimization, 2(4):649-664, 1992.
[34] Saverio Salzo and Silvia Villa. Inexact and accelerated proximal point algorithms. Journal of Convex analysis, 19(4): 167- 2012.
[35] Ruichen Jiang, Qiujiang Jin, and Aryan Mokhtari. Online learning guided curvature approxi- mation: A quasi-newtor method with global non-asymptotic superlinear convergence. arXiv preprint arXiv:2302.0885. 08588, 2023.
[36] D. Newman. Location of the maximum on unimodal surfaces. Journal of the ACM, 12(3):395- 398, 1965.
12

---

[37] A. Levin. On an algorithm for the minimization of convex functions. Soviet Mathematics Doklady, 160:1244-1222 1965.
[38] David B Yudin and Arkadii S Nemirovski. Evaluation of the information complexity of mathematicalpoogaammigg problems. Ekonomika i Matematicheskie Metody, 12:128-143 1976.
[39] Naum Z Shor. Cut-off method with space extension in convex programming problems. Cyber- netics and systems analysis, 13(1):94-995 1977.
[40] Leonid G Khachiyan. Polynomial algorithms in linear programming. USSR Computational- Mathematics and Mathematical Physics, 20(1):53-72, 1980.
[41] Pravin M Vaidya. A new algorithm for minimizing convex functions over convex sets. Mathe- matical programming, 73(3):291-333 -3341,, 1996.
[42] David S Atkinson and Pravin M Vaidya. A cutting plane algorithm for convex programming that uses analytic centers. Mathematical programming, 69(1-3):1 140, 1995.
[43] Yin Tat Lee, Aaron Sidford, and Sam Chiu-wai Wong. A faster cutting plane method and its implications for combinatorial and convex optimization. In 2015 IEEE 56th Annual Symposium on Foundations of Computer Science, pages 1049-1065. IEEE, 2015.
[44] Haotian Jiang, Yin Tat Lee, Zhao Song, and Sam Chiu-wai Wong. An improved cutting plane method for convex optimization, convex convxxe games, and its applications. In Proceedings of the 52nd Annual ACM SIGACT Symposium on Theory of Computing, pages 944-953, 2020.
[45] Bo Jiang. Haoyue Wang, and Shuzhong Zhang. An optimal high-order tensor method for convex optimization. Mathematics of Operations Research,
1000::119--11111 2021.
[46] Sébastien Bubeck, Qijia Jiang, Yin Tat Lee, Yuanzhi Li, and Aaron Sidford. Near-optimal method for highly smooth convex optimization. In Conference on Learning Theory, pages 492-507. PMLR, 2019.
[47] Alexander Gasnikov, Pavel Dvurechensky, Eduard Gorbunov, Evgeniya Vorontsova, Daniil Selikhanovych, and César A Uribe. Optimal tensor methods in smooth convex and uniformly convexoptimization. In Conference on Learning Theory, pages 1374-1391. PMLR, 2019.
[48] Yousef Saad. Iterative Methods for Sparse Linear Systems. Society for Industrial and Applied Mathematics, second edition, 2003.
[49] Martin Zinkevich. Online convex programming and generalized infinitesimal gradient ascent. In Proceedings of the 20th international conference machine learning (icml-03), pages 928-936, 2003.
[50] Zakaria Mhammedi. Efficient projection-free online convex optimization with membership ora- cle. In Proceedings of Thirty Fifth Conference on Learning Theory, volume 178 of Proceedings of Machine Learning Research, pages 5314-5390. PMLR, 2022.
[51] Kuczyáski and H. Wožniakowski. Estimating the Largest Eigenvalue by the Power and Lanezos Algorithms with Random Start. SIAM Journal on Matrix Analysis and Applications, 13444109941 1992.
[52] Aryan Mokhtari, Shahin Shahrampour, Ali Jadbabaie, and Alejandro Ribeiro. Online optimiza- tion in dynamic environments: Improved regret rates for strongly convex problems. In 2016 IEEE 55th Conference on Decision and Control (CDC), pages 7195-7201. IEEE, 2016.
[53] Peng Zhao, Yu-Jie Zhang, Lijun Zhang, and Zhi-Hua Zhou. Dynamic regret of convex and smooth functions. Advances in Neural Information Processing Systems. 13:12510-12520, 2020.
[54] Amir Beck. First-Order Methods in Optimization. Society for Industrial and Applied Mathe- matics, 2017.
[55] Alexandre d' Aspremont, Damien Scieur, and Adrien Taylor. Acceleration methods. Foundations and Trends@ in Optimization, 5(1-2):1-245 2021.
13

---

[56] Arkadi Nemirovski. Information-oen based complexity of convex programming. Lecture notes, 834, 1995.
[57] Yousef Saad. Numerical methods for large eigenvalue problems: revised edition. SIAM, 2011.
[58] Alp Yurtsever, Joel Tropp, Olivier Fercoq, Madeleine Udell, and Volkan Cevher. Scalable semidefinite programming. SIAM Journal Mathematics of Data Science, 3(1):171-200, 2021.
14

---

Appendix
A Optimal Monteiro-Svaiter Acceleration Framework
In this section, we present some general results that hold for the optimal MS Acceleration framework. In particular, in the first part of this section (Section A.T. present the proof of Proposition T In the second part (Section A.2).
we further provide some useful additional lemmas.
A.1 Proof of Proposition
To begin with, establish potential function for Algorithm | shown in Proposition 2 The result is similar to Proposition in [32), but for completeness we present its proof loosely following the strategy in [55 Theorem 5.3]. To simplify the notations, we use f denote the optimal f(x*)
Proposition 2. Consider the iterates generated by Algorithm 1 If f is convex, then
Ax+1(((Xx+2) - f*) + 2 p1llzt11 x²1 < (x( - fr) + 2 2/ - x x²1
(20)
Moreover, let a = 01 + a2 and we have
1 0² 1|z0 - x²1
X x yxll
(21)
Proof. Since f is convex, holds that
f(xx) 1(8++1) (V - 2 0, f(x*) - f(xx+)) - (V xx+1) 2 O.
By summing up the two inequalities with weights az and A₁ respectively, we get
A (K++1 xxx) 2 0.
Ax(f(x)) f*) (Ax + (f)(((()++))) - - (V 8c+1
ag
(22)
Let žx+1 = 8x+1 + Ak (x&+1 - xx). By rearranging the terms,
22 can be rewritten as
an
(A4 + ax) (f(3x1)) f*) An(f - f²) I (xx 26+1 x").
(23)
Moreover, note that the update rule for z in both (8) and ( can be written as
Zk+1 2₁ =
aA
V (x
(24)
7
Also, since also have Zς
Yk + As
(yr: xxx
from 2
write
as
Ak (8x+11 x5s)
Ak (yn xx:)
žb+1 Zb =
Xt+1 +
Yk
+
ax
as
A& + ak (fe+11 - yk))
ak 7i
(R&+1
yn),
(25)
where we used the fact that (Ag + (ax)m) M = am in the last equality (cf. 2b. Hence, combining (44 and 25 leads to
2x+1 Zktill = || Z/: (2k+1 Z1
as Tl
| - Vf mll S T |
4k
-y&ll-
(26)
where we used 7 in the last inequality. In the following, we distinguish two cases depending on nn 1/ or Ωp < Tik- In both cases, we shall prove that
1 A++1(f(XH+) -f) [z441 1=411xx111" < Ax(((xx) -f')+ 11z xx" -=xx-x11""
(1
oz)a
2n
(27)
15

---

If this is true, then Proposition 2 immediately follows. Indeed, since o < 1, the last term in the right- hand side of (27 is negative, which implies 20 Moreover, (2T) follows from summing the inequality in 27 from h = to N 1.
Case I: ûik = 7h. Since by we have XK+1 = 8k+1 and Ak+1 = As + ak.
23
becomes
( (x - (*) (x(s(x()) - (")
S ax (V((+)))))) 2 x*)
Using Zk+1
2/ aA V o (x in (8. we have
Ak 1 (f (x(x++1) f*) - (x( - f²)
Zk 1x+1 2x+1 x")
(z 2k+1
hx+1 xx+1} 1) + (Zx: 1k+1, Zx+11 - x²)
2
1kk 2441
11.
1.3
| 2k+1/1²
| 22 11?
(28)
1
7
+
Izz - -x"13
|izzz - -x"112
11--..11
2
2
1
(1
0²)
12
x
|38+1
x"1
||x++1
2
2
where we used (25) and (26 in the last inequality. This immediately leads (27 after rearranging the terms.
Case II: ii < Me- Since 0 < 7h < and Xx-1 = (1 (1-
ALtrees xx + 76(Aptan) & according to ( At Ax+7
by Jensen's 's inequality we have (As + (xak)[(xx+1) (Ax + 014 ([ + (1 7x.) A which further implies that
(As+ + ) (x(xx+1) -f (f(x)))
- f*).
Moreover, since Ak+1 As Jaak by ( together with (23 we obtain
A 1 (f(x f²) (x(x)x)) f*) S (V ēx+1 x")
Using 2k+1
2k TA4K V (((x+1)) in (9) we follow the same reasoning as in (88 to get:
A (f(x++1) f*) Ax(f(x)) - f*)
$ (z 1, 2x+1 - x*)
(2 2k+1, 1, 2x+1 1 1) + (7:: 2,+1, 2h+11 -x")
11 *x+1 11 1 + マ |xx+1 - -24+1ll
|| - 2/ 11
2
1
1
+ 2 |z x²112
|b x=1
2
(1 0²
|| - yx112, 2a)
S 2 1| - x
|6+1 - x1
211
which also leads to (27)
Next, prove a lower bound AN. Recall that B denotes the set of iteration indices where the line search scheme backtracks, i.e., B A {k: îp m...
Lemma 3. For any N 2 0, it holds that
An
N
Σ
1k
(29)
IEKEN
1,kgB
Proof. To begin with, according to the update rule of A++1 in and (9) and the expression of a in 2). the sequence } follows the dynamic:
AR+1 [ [Ax + agk if fip- = ?p (I: & B);
7k +
mh + 4m
where 7
and as
A + Yakk, if nt < (k E B),
2
Since we initialize Ao = O. we have dg 1/0- We further have A
1o. since we get A
Ao + dg no if 0 & B, while get A = Ao + 70a0
In no
no no if E B. Moreover:
16

---

In Case I where le $ B, we have
Ak+1 As + A& The +
1/
2 + Ab
(VAr + V Iik
2 A4 +
+
2
which further implies that V Ak+1
Ak + Likk
VAx Vis
In Case II where &s B, we have As+1
A& mkak 2 As. which implies that /A++1 V
Considering the above, obtain VAN 2 VA EISKEN-LAEB
which leads to 29
Lemma 3 provides a lower bound on AN in terms of the step sizes în in those iterations where the line search scheme does not backtrack, i.e., B. The following lemma shows how we can further prove lower bound in terms of all the step sizes {īis) k==
Lemma 4.
We have
Σ
V
1k
S
fiu
Σ
fK
(30)
v3
ISKEN-I,KER I,KEB
ISKSN -1,KgB
As a corollary, we have
/3 N-1
10
Σ
ik
(31)
2 /B
ISKSN-1,
1,4BB
k=0
Proof. When the line search scheme backtracks, i... k E B, we have iE S 37.. Therefore,
Σ
7h
Σ
BTG
BTp BBm=
Bm++1
(32)
ISKEN-
ISKEN-1KEB
we have 7 / B if k: B (cf. Line and M GM otherwise (ef. Line 13). This implies that 1j1 S no/ 3 and we further have
Moreover, in the update of Algorithm [
N
Bm
8T1+1
BBm=
Σ
V
Σ
B
k=1
ISKEN-2KEB
ISKSN-2
2,KEB
S vmn
Σ
Σ
V
Biij
ISKE 2
1CKSN-2,KEB
< Viuu +
Σ
VÍN +
Σ
VBin-
(33)
1ckEN-1,k$B
ISKSN-1,KBB
We combine
32 and 33 to get
Σ
Vin
V
Σ
7JK
Σ
BBn
ISKSN J,K
ISKSN
I
By rearranging the terms and simple algebraic manipulation, obtain (30) desired. Finally, 31 follows by adding vij + E Vîj. to both sides of (30
Now are ready to prove Proposition
Proof of Proposition y
By Proposition the potential function øu &
|
non-increasing in each iteration. Hence, a recursive augment have A
Φv < Φu l|z 12 which yields f(xx) < || -x112 Moreover, combining
24
Lemma 3 and 31 in Lemma leads to the second incquality in Proposition
17

---

A.2 Additional Supporting Lemmas
crucial part of our analysis is bound the path length of the sequence [y This is done in Lemmaß8 To achieve this goal we first present the results in Lemmas 5 which provide the required ingredients for proving the claim in Lemma In our first intermediate result, we establish uniform upper bounds for the error terms
I x"1l and ||x x'|
Lemma 5. Recall that 01 + 02- For all & 0, we have |
x"1l I 2|zo
x'1l and
|xx: x
Proof. To begin with, it follows from 20 in Proposition that
1 2 1z0 112
2 Iizx 11 - -
2k x 112 S ff(
=
Hence, we get 1zc - 1| S |z0 || for any & 2 0. To show the second inequality, we distinguish
two cases and in both cases we will prove that
az As)
Ak+1 [x11 x11 1| < x|
| yz 1x 116-44+11+61111 z6+1 1x=1 x=1 As]1|z++1
Mh
(34)
Case I: flk = T]k- Recall that in the proof of Proposition we defined 2x+1 = * + As au
(fx+1 - xxk)
Since X++1 X6+1, we have Xt++1
Aatax
As
Xb + Acto,
Zx+1 and by Jensen's 's inequality
| - x111 S A
||x - x" 112
|ī - x x"11²
A +
A + a
Furthermore, we have
20²a m
| x*1 S 2||2++1 Zk+1
1 yx -x'13, (35)
where we used 26 the last inequality. By combining the above two inequalities, we obtain
112 < Axl|xx x" 112 + 20
(Ac + <x))| X
[ y&ll2 + 2ax [ - -x113
which leads
34 (note that A A aa in Case I).
Case II: Since X
Al- Autibaa
Xx +
(((++oo)) Sx+1 and 8x+1 Actias
Ak X Actas ăc+1, ANtay
have
A
7akk 2o+11.
X
Xx +
A makk
An + 7B0k
Similarly, by Jensen's 's incquality we have
(Ax + 7x0+)| - x11 < A - -x111= 110111221111 - -x1112
Combining this inequality with (35) obtain
20 a
(Ag+77aa (a(ak)))))))) |[xx+1 1|" S A& X x"1 1|2 a.
(36)
mh
which leads to
34 (note that Ak+1 A& na in Case IT).
Now by summing over h: O,..., 1, we get
N-1
20
N-1
An Xx x <
(AR+1 - Av)
| 1
2(A+1 Aw]lzz+1 - xx112
nk
0=
k=
(37)
N-1
N-1
| - yxll 2||zo xx x11 N-1
< 20² 2
(AR+1 As)
Σ
Σ
(Ap+1
AA))
m
(38)
20²
AN ||z - x=1 2.A ||z - 11
(39)
1- σ²
2Ay
120 x²1
(40)
1 σ2
18

---

Hence, this implies that ||x: x² x*11 11 S 1=zz - x²11 for any k 2 0.
key term appearing in several of our bounds is t In the next lemma, we establish an upper bound for this ratio based on a factor of its previous value, for both cases of our algorithm.
Lemma 6. Without loss of generality assume /3 1/5. Case we have AusitRess Otherwise, in Case II we have
v3 Ag +ax
1611 S Ak41-0k11
2 at=.
vB+1 Actor"
Proof. By the choice of ag in we have (Ax + axk) = az for all to 2 0. As a result, we have
0
27
2
A +
a%
T] v 1 + An
1 1+
4
A
and similarly
ak+1
2
Ak- ax+1
1+ + 4
Act1
7
In Case I, we have m 7 and Ak+1 2 Ag. Hence, it implies that AR+1/16+1 2 3AA/ which leads to
0x+1
1
2
a2
A&+1 + ax++1
+
48A
vB B +
48Ar 7
V3
1+
1
4 n
vB A + ap
where the second inequality follows from the fact that B S 1.
In Case II, we have M7++1 = ii
7A77 and ARt1 = AA + Okk Since we also have ag 2 7) and Th S B, we obtain Ak+1/ / (mun)) +12 A&/ A + Hence,
a&+1
2
2/B
2
2/B
<
Ak+1 + 0g+1 1 V
55+ +
1Aa Di
1
1A
vB 1 +
V
1AE
vB +
A ax
VD
where we used B > 1/5 in the second inequality and the fact that 2 in the last inequality.
v3" x 2
10+1 2/3
(1 + x) for
Remark 5. 3 - 1/5, then in Case we still have Ast1 I 1 at
while in Case II we a < V Actis Thus, in the case where B S 1/5, the derivation below still holds have Ax11+a except that the absolute constant C2 will be different.
ax+1
vB A,=a;"
og
Next, as a corollary of Lemma we establish an upper bound on the series N-1 - Moreover, we use this result to establish an upper bound for N ||*++1 - yull-
Autar
Lemma 7, We have
AN 1 2/3 3 - log A )
a
(41)
A& + ak
vB 3 3
A=
Moreover.
1 1+ 8
Av
| - S
V o2 -3-3
1 +
log
A
1:20 xxIl-
(42)
k=0
Proof. Given the initial values of A₁ and og we have
1
aj
Σ
Ok
(43)
A + ak
A + as
A as
Au + 0s
k=1
kEB
k
19

---

Note that using the result in Lemma
6
N
Σ REBAA21
Ok
0k
(44)
Ak + ak
A
+ ax+1
Σ
a&++
E
0K+1
(45)
AR+1 ax+1
Ak+1 + ak+1
kEB
kEB,120
Σ
0 V3
Σ
2/B
Qk
(46)
A
ak
VB +
A2 + O2
kEB
kEB,k20
Ok
2 3
ag
(47)
VB
ak
VB 1 A ax
k4B
Hence, if we move the last term in the above upper bound to the left hand side and rescale both sides of the resulted inequality we obtain
1
3
1+
+ a2
VB
A + a
KEB
kgB,A21
Now, if replace the above upper bound into 43
obtain
N
1 2/ 3 1+
Σ
Σ
(48)
A + ak
vB
3
A + ap.
k=
REBA21
Moreover, note that for lo & B, have Ak+1
A Or. Hence,
Gg
Σ
As
(log(A++1)- - log(An))
A +
Ax+1
R
R
N-1 1 1
(log(Ar+)) - log((A)) = log An
A
Now if we replace the above upper bound, i.c., log Ax with EeeB, Axtax into the expression in the right-hand side of (48) we obtain the result in (41). Next, note that by Cauchy-Sehwarz
Schwarz inequality, we have
N
N-1 N-1
N
-
- y&ll S
[
yal2
Zo
x"11
R==
k=U
where the last inequality follows from (2T). Moreover, based on the expression for in and the result in 41
that we just proved, we have
N-1
1 + 1+2 - 336 3
An
Σ
vB - 3 1+ log
(A
a
A + ak
A
k=0
k=0
k=0
Combining the two inequalities above leads to
(42)
Now ready to present and prove Lemma 8 which characterizes bound the path length of the sequence {yx) Lemma 8. Consider the iterates generated by Algorithm 7
Then for any N.
N-1
AN
Iyy+t1 yxll S C2 1
log
1220 -
xx x-11
A
K=0
where
1 22
1+B - 3
1
2
1 + 2/3 -
C2 = 2/
V3 -
+
VB
1 +
V
1
VB B
(49)
20

---

Proof. By the triangle inequality, we have
Vy yk+ill & |8++1 - yxll ||&t+1 yk+ill-
(50)
We again distinguish two cases.
Case I: iis- The In this case 88+1 X1+1 and yk+1
Aut1 Xx+11
ac+1
2k+1 hence
A4+1 A ax+1
|| yk+1 = ||x+1 -k++il a++1llza+1 Xx+1 ||
1
0x| Izo - x²
1+
Ak+1 ak++
1 x²
Ag + ak
where we used Lemma and the fact that || - X
|a+1
x" +
|X++1 x S
(1+ + V )11o - x=1 I| in the last inequality. Therefore, using (50) and the above bound we have
I1v yktill s |x++1 yxll + vB
1 +
2
0.
Iz - xx1
(51)
1
1 0²
A& ak
Case II: ilp
Th-- Since Xx+1 Autrua
Av
X +
Ttap Aut
%x+1 and Xx+1
Autay
Xx +
ap Aotas
2b+1-
we get
As A + ak
(x
Talk (xx+1 (x*+1-2++1))+
Ok
Ax + 7
x++1+ (1 7x)aa 2k+1-
S6+1
2
Ak
A + Ok
AB + ax
Ak + ak
Thus, given the above equality and the expression Yk+1 have
A212+aR11 Xk+1 + Axi1tax Zk-+1,
Ah11
aNll
| - yk+1ll S
(1 (1-7ks)ak
||
Zk+ill +
(1 - 117k)a
0k+1
2k411 Xx++1 ||-
A& ak
A + ak
A&+1 + ak+1
(52) yE
Moreover, based on the result in 26 which implies that
upper bound ||ăx+1 2x+1 || by
th 7a
||
(1 - 7 (1 - -1%kaa [ Zu+l S 0
|| yxlll = o(1 - (x) 1x yvll S || Vkll
Aş ag
1)s (A4 ax)
where the equality holds due to the definition of ak, and the last inequality holds as both y and in (0, 1). On the other hand, note that
)aa
0K+1
(1 7R)aR
(53)
A + ak
A + ak+1
A + ak
A + ak
0g+1
(1 7k )k)ak
2/Bax
(1-7k)a
ag
<
(54)
A + 01+1
A + as
vB + I(Ar + +ax)
A + as
A ak
where in the second bound we used the result Lemma and the fact that 2 vD+I
< Hence, we
get
| yk+1 < || ak llzx411111111 Xa+Ill S |8411 1-841-yyll+++ 1+1
0 ax||z0 |20 x"1
A + dx
o²
Aş + ag
where the last inequality follows from the fact || x S the bounds Lemma Now by applying the above upper bound into (50)
1k+1 x" || + | - x" || and
we obtain that
llyx yk+1ll < 2 2
yxll + (1+
2
ak
Zz
x ||-
(55)
o2
A + ax
Considering the upper bounds established for lly& Yk+1ll || in case (equation (equation (55 we can conclude that
51 and case II
1
a2 A + ag
Iy Yetill |I $ 2 2|| yxll
1
|
x+|
(56)
V3
1 σ²
Finally, Lemma 8 follows from summing (56) over lo = 0 to N - 1 and the result of Lemma 7
21

---

Subroutine Backtracking line search
1: Input: iterate y I2 gradient g R.. Hessian approximation B E S initial trial step size v)
2: Parameters: line search parameters / E (0, 1), 01 2 and 02 0 such that Q1 02
3: Set f) f-1 ?
4: Compute S. + LinearSolver(1 + AB, -ig; an) and &- + +
5: while |x+ - / ||2 > (ar + az)ll$+ y112 do
6: Set
x4 and 3n
7: Compute 5t LinearSolver(I rjB, -B; an) and x f- y 8
8: end while
9:
if then
10: Return 1) and X+
11: else
12: Return ij. 8+ and x
13: end if
B Line Search Subroutine
In this section, we provide further details on our line search subroutine in Section 3.1 For complete- ness, the pseudocode of our line search scheme shown in Subroutine Section B.1 that Subrountine will always terminate in finite number of steps. In Section B.2 we provide the
we prove
proof Lemma
B.1 The Line Search Subroutine Terminates Properly
Recall that in our line search scheme, we keep decreasing the step size v by factor of 3 until we find a pair x4) satisfying (IT) (also see Lines and in Subroutine []. In the following lemma, we show that when the step size is smaller than certain threshold, then the pair (f), x+) satisfies
both conditions in (T0p and (1 which further implies that Subroutine I will stop in finite number
of steps.
Lemma 9. Suppose Assumption 7 holds. if < +1+ op
and & is computed according to 12
then the pair (î) x+) satisfies the conditions in (10) and 111
Proof. By Definition] I the pair (f) x4) always satisfies the condition in (1O) when x is computed from (2).. Hence, in the following we only need to prove that the condition in (TT) also holds. Recall that g Vy()). By Assumption I the function is L1-smooth and thus we have
[V gll = [[V ff(y)) s L1|| yll.
Moreover, by using the triangle inequality, we get
|| f*++) B( - y)ll < [V - gll + B - y) s (L1 +
[Bb+p)] yll.
Hence, if f S T. Bop
have
fx++) - g - B(* - y)1l < 02 |8 - yll-
(57)
Finally, by using the triangle inequality,
combine
10 and (57 to show that
| yy + fx++))
| y + f(g + ( - y)) + fx*+)
B( y))ll
ب y fig B( y)) +|0||+10((1/ (I B( y))l
< 0n|/5 - y +2 yll (al +a2) -
which means the condition in 11 is satisfied. The proof is now complete.
B.2 Proof of Lemma
We follow similar proof strategy as Lemma 3 in 35 In the first case where A $ B, by definition, the line search subroutine accepts the initial step size T)k i.e., nn = mh- the second case where & E B. the line search subroutine backtracks and returns the auxiliary iterate 2 which is computed from
22

---

12 using the step size îix - ft./B. Since the step size ÿt is rejected in our line search subroutine, it i.e.
implies that the pair (X++1. Ma)) does not satisfy
[ Yk + 7 (x
(a1 a2) 8x+11 ykel--
(58)
Moreover, since compute X from 12 using step size tjs the pair (īn, & also satisfies the condition in 10P. which means
| Yk + im f(yx) + B& (x - ya))ll - a. || yull-
(59)
Hence, by using the triangle inequality, we combine 58 and (59) to get
V vf(y))) (x++++1 yx)ll
| yk + 7/ (x**++1))++++ |x Ya The (Vf(y)) + Bx - yk))
(a1 + a2) 1 yxll || - 01/ 0/|/8+11 - -vall
02 Sx+1 yull,
which implies that
a2 02 - -vull
Bîj
f Vf(yk) - B2 yk)ll
This proves the first inequality
13
To show the second inequality in 13 first note that X and SR+1 are the inexact solutions of the linear system of equations
îik Bx)(x yk) -j and (I Bx)( - yk) -
respectively. Let x and 84+11 be the exact solutions of the above linear systems, that is, x 1 (I
y: nu (I + în Bx)
3x))-ge gs and x Yu
ijs Bx) -1ge- We first establish the following
inequality between || yull and |8 vx
1 3
|| yxll
| -kll--
(60)
This follows from
[ y& = | (I+ gk $ II+ Ba) -gxl
|
ykll,
where we used the fact that (I + ij B+) 11 = (I + 1) Bb) -1 in the first inequality. Furthermore, we show that
(1 01 a1) [ y < | yıl < (1 + [n) yull,
(61)
(1 01 [1]]| 8 - y < | yxll < (1 + 01) X - Yall-
(62)
We will only prove (6T) in the following, as satisfies the condition in (10). we can write
62 can be proved similarly. Note that since (fx, 8x+1)
| Yk + (g+ + BB - y&)) = [I + Bx) - x S a1 [1 yxll
Moreover, since B& 2 O, we have || - x S |I ( x < 01 | - yxll-
Thus, by the triangle inequality, we obtain
[8 y&ll < | Y&ll
| - S || S (1 an)) ykll-
|| || - yxll 2 || yvll
| - S || 2 (1 a an) X++1 -vll--
which proves 6II. Finally, by combining (60). 61 and 62 we conclude that
|x+1 yall < 1
1 1 (1 a1)
| yxll
x yxll -
| -vlll-
α1
(1 a1)
This completes the proof.
23

---

C Hessian Approximation Update
In this section, we first prove Lemma
in Section C.1 and remark on the computational cost of
Euclidean projection in Section C.2
Then we present general online learning algorithm using
approximate separation oracle in Section
C,3 and fully describe our Hessian approximation update in
Section
C.1 Proof of Lemma
We decompose the sum
N-11 kk=0
as
%
N-1
Σ
1
(63)
22
k=0
1
ISKEN-1,KEB
"
1SkgN-1,kdd
Recall that we have im = 7) for & B. Hence, we can further hound the last term by
Σ
1
Σ
4K
IEKEN-
EEEN
1,8
k=1
1
Σ
Σ
+
7
Mh+1
ISKEN
2,kEB
1SK$N-2,kg
Recall that we have 7 = nK if le E B and m6+1 = n//3 otherwise. Hence, we further have
Σ
Σ
1
Σ
+
1/1
1h+1
2
ISKEN-
1,
ISKS
KEB
ISKSN
k,
B2
Σ
|/|
f
1%
I<KEN
ISKEN -2
B
Σ
7
1%
ISKEN-
TEREN
1
By moving the last term to the left-hand side and dividing both sides by - B2, we obtain
Σ
I 1
(64)
32
qi
ISKEN 1
ISKEN
I,KEB
Furthermore, since ?/1 2 no. have m
- Hence, by combining
63 and (64
we get
N-1
2 82
1
+
Σ
1 )
2 - 12 (1- 32)0
+
2 g-
(65)
1 82
1h
1- 32 OSkSN- 1,KB %
0=
where in the last inequality we used the fact that îj oo if 0 & B. Finally, (1)) follows from Lemma
and
C.2 The Computational Cost of Euclidean Projection
Recall that Z P B E Sd 0 = B = LJ}. As described in [35 Section D. the Euclidean projection on has a closed form solution. Specifically, Given the input A E g we first need to perform the eigendecomposition A VAVT. where V is orthogonal matrix and A diag(A),...:AA)) is a diagonal matrix. Then the Euclidean projection of A onto Z is given by VÁVT where Á is a diagonal matrix with the diagonals being Å = min(L,, max (0, )x)} for 1 S 1: < Since the eigendecomposition requires ((d3) arithmetic operations in general, the cost of computing the Euclidean projection can be prohibitive.
24

---

Algorithm 2 Projection-Free Online Learning
1: Input: Initial point wo E Bn(0), step size > 0, 6
2: lor(= 0..... T
do SEP(W1; 6z)
3: Query the oracle (Yr,
r I then Case we have w, E
5: Set t Wz and play the action Xt
6:
Receive the loss (xx) and the gradient g: = Vxx((x))
7:
Set 8:g:
8:
else Case H: have wa/"t E C
9:
Set X f- we/7) <We///
and play the action Xt.
10:
Receive the loss (((x;) and the gradient g: = (x,(x.)
11: Set g. + ge + max [0, - (g). xs)
12:
end if
13:
Update W maxi w - PRI R3 (Wi pg.)
Euclidean projection onto Bg(0)
14: end for
C.3 Online Learning with an Approximate Separation Oracle
To set the stage for our Hessian approximation matrix update, we first describe a projection online learning algorithm in a general setup. Specifically, the online learning protocol is as follows: For rounds observes a loss function & R"" + R. We measure the performance of an online learning algorithm
.. T o....... 1, a learner chooses an action x; E C from a convex set C and then
by the dynamic regret 49 122 defined by
T-1 UT Σ (x,) Σ E(le))
D-Regy (u1,
where {ur}+= a sequence of comparators. Moreover, we assume that the convex set is contained in the Euclidean ball (0) for some R > 0, and we assume 0 E C without loss of generality.
Most existing online learning algorithms are projection based, that is, they require computing the
Euclidean projection on the action set C. However, as we have seen in Section C.2 computing the
projection is computationally costly in setting. Inspired by the work in (0)
will describe an
online learning algorithm that relies on an approximate separation oracle defined in Definition 3
Definition
The orucle SEP(w; 5) takes w E Bn(() and > as input and returns scalar r > and a vector E Thn with one of the following possible outcomes:
Case I: y I 1 which implies that w E C;
Case II: 7 > 1 which implies that w //1 E C and (s.w-x)27- 27
Vx EC
To sum up, the oracle SEP( (w; 8) has two possible outcomes: it either certifies that w is feasible, i.e., w E C, or it produces a scaled version of w that is in C and gives an approximate separating hyperplane between and the set C.
The full algorithm is shown in Algorithm 2 The key idea here is to introduce surrogate loss functions (w) = (gr, w) on the larger set (0) for S l < T 1, where g: the surrogate gradient to be defined later. On high level, we will run online projected gradient descent with (w() to update the auxiliary iterates (we)er heo0 (note that the projection (0) is easy compute), and then produce the actions {x, Fizu for the original problem by calling the SEP (we; 5) oracle in Definition 3 The follow lemma shows that the immediate regret (e(w;) - (x) can serve as an upper bound on (x,) (x) for any x C.
Lemma 10. Let {-
1T-1 be the iterates generated by Algorithmm 2 t t= 0, 1, 1,..., T 1. Also, for any C. we have
Then we have X E for
(gt, X4 x) v (gc Wi x) + max(0, (gg, x,) 15,
(66)
1
|w x1²
We+1
x112 + S/2 1Brr12 max(0, (0,00080 (ge, xx) J8,,
(67)
2p
2p
and
Lg S |g: (g. |s|1|
(68)
25

---

Subroutine Online Learning Guided Hessian Approximation Update
1: Input: Initial matrix Ba S s.t. 0 I Bo I L, step size p 0. 0. (+))
2: Initialize: set Wo + Ti (Bo 211). Go + L1 (B) and Go +- Go
2
3: for 1,..., T
4: Query the oracle (Tr. St) < SEP(WW;;,,)))
5: ie S 1 then # Case
6: Set B: + w, and B, + Ln B: + I
2
7:
Set G. + 2, (B.) and G. + G:
8:
else
Case
9:
Set B: + Wa/7 and B. + L B. + 1
10: Set G. + 2, (B.) and G. t G. + max (G,, B.))S,
11: end if
12: Update W441 max[va.) GGlle) (W.-pG.)
Euclidean projection onto B, Baa(0)
13:
end for
Proof. By the definition of SEP in Definition 3 we can see that x, E C for all 1...., T.. We now show that both (66) and (68) hold. We distinguish two cases depending the outcomes of SEP(w,; 6.).
* If ya S 1, then we have = We and g: In this case, (66) and (8 trivially hold.
If 7o > 1, then Xx = wt//t and g: = gt + max (0. - (g:, x1) We can then write
(gr, Wt - x) (g: max(0. (g. x /s, W: x)
Yx x) + max 0. - (gi, xx) (s,, W. - x)
(gr, X2 x) + (* (g., + max - - x,) (7) - 6,) (gr,
X; - x] max (0, (g. x, 16g + (% 1) max(0, (g.. x1) (g., X x) max(0. (g., x1) 36b.
which leads to (66) after rearranging. Also, by the triangle inequality obtain
g: < ||g:|| + ma - (g. xx) | < ||g:l| + (g)
which proves 668
Finally, from the update rule of We+1, for any x E C C (0( we have (w, - pg: -W+1, W x) 2 0, which further implies that
(gr, we x) < (Br, W. W:+1) + (W: - We+1, Wo-1 - x)
(69)
1 ||, w 112 (70)
(B:, W. W
|| -x122
|/111
x122
2p
2p
2p
W: x|
|w++1 - xx12
+ p 2 B: 1.
(71)
2p
2p
Combining (66) and
(71) leads to (67)
C.4 Projection-free Hessian Approximation Update
Now we are ready to describe our Hessian approximation matrix update, which is specific instantiation of the general projection-free online learning algorithm shown in Algorithm 2 In particular, we only need specify the convex set well the SEP oracle.
Note that we have = [B E S 0 < B = L in our online learning problem in Section 3.2 Since the projection- scheme in Subroutine Sartuinnerrrqee requires the set C to contain the origin, we consider the transform B & 2 [B-
2 I) and let C 2 e [ Sd -I B I)- = (B e sd | lop <
1}. We note that 0 E Z and E C B,) (W edd || S vdy. Moreover, we can see that the approximate separation oracle SE; ā, g) defined in Definition D
to a stochastic version of the oracle in Definition The full algorithm described in Subroutine 2 and we defer the specific implementation details of SEP(W;& 6. q) to Section B.2
26

---

D Proof of Theorem 1
Regarding the choices of the syper-parameters, we consider Algorithm| with the line search scheme in Subroutine [ where a1, 02 E (0. with + 02 < 1 and 3 E (0, 1), and with the Hessian 2 where = 128
1/2,5(++1) log (++1) for t 2 1, and 2)) for t 2 O. In the following, we first provide proof sketch of Theorem The complete proofs of the lemmas shown below will be provided in the subsequent sections.
approximation update in Subroutine
6; =
Proof Sketch. To begin with, throughout the proof, we assume that every call of the SEP oracle in Definition 2 successful during the execution of Algorithm Indeed, by using the union bound, can bound the failure probability by Σ lt S P
[
2.5 Et=2 Elog" S p. In particular, we note that
Subroutine ensures that < BB 3 I for any k 2 0. 2
We first prove Part (a) of Theorem I which relies on the following lemma.
Lemma 11. For B, we have 4 (Bx)
We Bu
$ LT
15
We combine Lemma and Lemma 11 to derive
3- +
2 2-3 -2
Wx BySK 1|?
2-
(2 (
+
N.
(1- 32)0
(1 (1-33 (2)33))
[sp 12
-1- (1-32)or 32)0g
(1 3)
REB
By further using (14) and the elementary inequality that Va Va + vb. we obtain
CaIa
|20
1²
Czl1zo x²1
f(xN) f(x')
N
No.5
where C4 = CIV (1-8-jøg
2-32
2- (I- PI-B3aa32 and C3 Ci
(2- (2-82)
11-3"
Next, divide the proof of Part (b) of Theorem 1 into the following steps.
We first regret analysis to control the cumulative loss ΣTs (B,) incurred by online learning algorithm in Subroutine 2 In particular, we prove a dynamic regret bound, where we compare the cumulative loss of our algorithm against the one achieved by the sequence {H.)
Step 1:
Lemma 12.
We have
T-1 & (H,) + 2 8 + 512L
T-1
T-1
(B,) S 256 B.
Holls +
Σ 1=0
Vd
H,lr,
t=0
1=0
where H. (yy.)
Step 2: light of Lemma
12 suffices to upper bound the cumulative loss ΣT- & (H,) and the
path-length Σ(=1 |H
H,l in the following lemma. To achieve this, we use the stability
properties our algorithm
21 and Lemma which most technical part of the proof.
Lemma 13. We have
T-1
C3 4
AN
Σ
((H)) <
121100--
x
and H Helle S C2VdL2 1+log
||z 1x0-x"11 x*11.
A1
(72)
where C2 defined in
(1+a1)² 49 and C3 g2 (1-01)((1-02)
Step 3: Thus, we obtain upper bound on ET- ((,)) by combining Lemma 12 and Lemma
Finally, in the following lemma, prove upper bound on I by further using Lemma and Proposition] I
Ag
27

---

Lemma 14. We have
max Li
and ou 1 v VM
M L2dd x log
AN
N2.5
where we define log" (ar) max (log(z))))
Ce
M C,L + CB Holli CgL x=11 CioLi CioLi x"11
and C₁ (i ..... , 10) are absolute constants given
Cs 42 1 B²) 32 C7
5C 033 32
256C 033
C3Co 0330
512C2CB 03.32
Cg
C9
C10
Therefore, Part (b) of Theorem immediately follows from Proposition
I
In the remaining this section, we present the proofs for the above lemmas that used prove the results in Theorem
D.1 Proof of Lemma
(((g++) V(yx)) and S& S Ye for k E B. We can write [((*++1) Vy(y))
Recall that Wk
($++1 - Y&) by using the fundamental theorem of calculus, where H& fo 2 f(xx++1 + (1 t)ys) dt. Since have 0 f2 3 L for all Rd by Assumption 1 it implies that 0 < Hg I. Moreover, since 0 N B L,I, further have -L,1 = He B N L,I, which yields ||H Ba I L1. Thus, we have
w BB = [| Bx) yk)ll S L 1x ykll.
which proves that lx (Ba) LZ.
D.2 Proof of Lemma 12
To prove Lemma 12 we first present the following lemma showing smooth property of the loss function 4.. The proof similar to (35 Lemma 15].
Lemma 15. For k B, we have
1 (-sx(wx B8²) (wy, Bsx)ss)
(B)
(73)
|s
Moreover, for any B E Sd, it holds that
[V) S VV (B),
(74)
where || ||e and
|| denote the Frobenius norm and the nuclear norm, respectively.
Proof. is straightforward verify the expression in 73 The first inequality in 74 follows from the fact that || S || for any matrix A E $d. For the second inequality, note that
& [I. S
| (I|sk (W2 Bsa) T|. + ll(w B8&)s (1.)
2|| Bsa
| [W
Bs BssI|||sh|]
2 (B),
||sk||
where in the first inequality used the triangle inequality, and in the second inequality used the fact that the rank-one matrix uv has only one nonzero singular value |
[all|v|l
We will also need the following helper lemma.
Lemma 16. the real number satisfies - A Bvz, then have S 2A + B2.
28

---

Proof. From the assumption, we have
B
B
VI
S
AA+
2
4
Hence, we obtain
2
1
B
x S
[V++
+
S 2A + B²
4
2
Before proving Lemma 12 we also present the following lemma that bounds the loss in each round.
Lemma 17. For any H E 2, we have
(B.) < 46, (H) 64L]/ H 6642 64L7 HI + +22,87
Proof. By letting x, B, x= H
2 H L).. G. ~ 2, (B,) g. -Go, Wt W.
in Lemma I0 we obtain:
(i) B, E 2, which means that [Blllop S 1.
(ii) It holds that
(G., B.
H) S 2p
W. H 2p
W - H + G.G + max -G,, B,) Jo,,
(75)
||&.|l= - |.,|| F + (Gr
(76)
First, note that || < 3 by Definition (76). get
and (G,, B,) S
S GGG] Together with
16 (B,),,
G |r - |G + 3/ < G. S
(77)
L1
where we used the fact that G, Lı since l, is convex, we have
2. (B,) and Lemma 15 in the last inequality. Furthermore,
L4 2
(BB.) - 4(H) < (V), H)
(G., B: H),
where used G; (B,)) B, & L combining 75 and
(B. 1). and H (H II).
Therefore, by
get
[
T
1 |. (78)
4(B.) - (H) -
W. H1
W (+1
H + LT G.
112
8p
8p
4
I
I
<
W.
HI
W
(+1
H| + 32pb, (B.)) L1 (B 8
(79)
8p
8p
Note that (B.) appears on both sides
77).
By further applying Lemma 16 we obtain
L
(BB) < (H)
W.+1
Á + 64pl, (B.) + L788
4p
4p
Since p 1/128, by rearranging and simplifying terms in the above inequality, obtain
& (B.) < (H) 6,
H
- H + 2L788.
29

---

Proof of Lemma 12 let H, = f(y) for
T .. 1. Thus, we get
((B)) I (H, + | W. H.1 64L
W L+1 H + 2L
= -H, + 647,] W. -H,1 - 647, W i+1 - -H++1li 2T,76
+ 6444 (1[ 1+1 H.+1 -
Furthermore, note that
We+1 H 112 | 1+1 H,1l?
-([W++1 He+1
+
W H:Il=) - Hetille [ -H.lF)
8vd
SAV0|IH H L
H4+1 H.lle,
where in the last inequality used the fact that H. H441, W4+1 Bva(0) and the triangle inequality. Therefore, we get
(B,) s (6, -.W 6444/ W 2+11 H. 1
H,IFF.
By summing the above inequality from to T we get
T-1 Σ (B,)) 64L7 Holl p +
T
T-1
T
& (H,) + 21,
Σ
6, + 5121
Vd
H - H.lle-
E=
1=
1=0
Finally, we use the fact that Wo A b
(Bu
L I). and Hy
A In
2
HH 1) I) to obtain Lemma 12
2
D.3 Proof of Lemma 13
By Assumption 2 we have |w, H,s;|] [V(((+1)
- vf(y)) Vf(r))
yo)l 5
La Re+1 yell.. Thus,
[w, H
L2 4
(1 a1 L23
4 (H))
8e+1 - y2
4 )2 |8411 yrll?
|s 1
where we used Lemma in the last inequality. Also, Since as 2 ni for all h: 2 0, by
21
we get
N
-
Zo x
Hence, we have
(1+a))21) 412 -a1)²
(1+a1) 2 N-1
HH)) <
Σ
| - yxl
Σ
|x++1 yx
4. a11) k=
AEB
(1 2 12 x11
45 - 21) 21 02)
which proves the first inequality in (72). Furthermore, by Assumption 2 we have
|| == f(y=l v² f(y:)l f /(/l=
- -vll.
Hence, by using the triangle inequality, we can bound
T-1 H. Helle S VdL2
N-1
Σ 1yk+i - yall VdLzcs (1 1+gog
Av 1
|20 x1
where we used Lemma in the last incquality.
30

---

D.4 Proof of Lemma 14
Before presenting the proof of Lemma we start with a helper lemma that shows lower bound on A1.
Lemma 18. We have A1 ijo 2 oo T. J.
Proof. The equality A1 no shown in the proof Lemma 3 To show the lower bound on nj. we use Lemma and separate two cases. If 0 & B, then we have no = co. Otherwise, if 0 E B, then we have
82 yoll
[ V
yoll
Moreover, as shown in the proof of Lemma 11 we have [ f(yo)) - Bx yooll -
L || - yoll. which further implies that fju- 2
oaf This completes the proof.
We combine Lemma 12 and Lemma 13 to get
Ww
B
(B,) s 256 B. - Holl + CaLă - -x"11 2L3
KEB
t=
AN
+512C2L, 1 log
|zo x'11.
A1
Since 8, = 1/(lv + 21 + 2)), we have
T-1
Σ
Σ
o
tn²
T+1
th2
1
2
tln t
dt +
2 2
In
1
In(T + 1)
2.5.
t=2
Hence, it further follows from [4
and Lemma 2 that
N
N-1
S 4C²
Σ
A
kh
k=
4C 332
4C 82)
W
B&SK
(1 P00
|s
REB
Co + C+Li + Os| -Hol? 2 Co LoL3IIzo x11
Av A
+ C1obi L2d
1+
)
1120 xx11.
(80)
To simplify the notation, define
+ +
Holl? Cg x21 112 + CioLi x"I1,
and the inequality in 80 becomes that achieves the equality
N 43
< M L1
x log Av A
Let A be the number
N M CioLy xx1
A log
(2
A
and we can see that Av AÑ. Thus, we instead try to construct lower bound on Ai. If A S A then log(A;;/A1) I 0 and furthermore
N
VM
)2 S S
(81)
Ax
N2.55
Otherwise, assume that A3 > A. Then log(Az/ AA)) > 0 and we first show an upper bound on Azy:
N5
A
= M + Lodol - x" log
2
M
AW <
22.5
)²
A
31

---

This in turn leads to lower bound on Ay:
N
AN
max
Li
JNN5
(A M + C1oLL 1l1
M + CIoL1Lad1z0 x'1 log
and do
Ai
VM
where we also used the fact that A 2 min azd) } (cf. Lemma 18
Thus, we get
max 1
J22.
1
M Cio Lzd|zz
x"1
log
ag og
(82)
Av AN
N2.5
M
Combining both cases in
81 and
82).
conclude the proof of Lemma 14
32

---

Subroutine 3 LinearSolver(A, (A, a)
Input:
$4 E R.
2: Initialize: So < 0, rp 6 Aso. Pa t- ro
3: for k = 0, L.... do
4: if ||rl|l22 S a||5s|12 then
5: Return 8k
end if
7: ar 4 (Fk, Art) (Ape,
8:
re+1 t- Fa a4 Ap;
10:
Compute and store Ar
11:
3e + (I4+1, Art-1) Are)
12:
P++1 P441 BLPK
13:
Compute and store Apk+1 Ars+1 + App
14: end for
E Characterizing the Computational Cost
In this section, we first specify the implementation details of the LinearSolver oracle in Definition) and the SEP oracle in Definition 2
Then in Section E.3 we present the proof of Theorem
E.1 Implementation of the LinearSolver Oracle
We implement the LinearSolver oracle by running the conjugate residual (CR) method 48 solve the linear system As = b. In particular, we initialize the CR method with Sp = 0 and returns the iterate S&
achieve || Ass bl| S a/|ss||-- The following lemma provides the convergence guarantee of the CR method, which will be later used in the proof of Theorem
Lemma 19 (56) Chapter 12.4)). Let s* be any optimal solution of As* = b and let {sk} be the iterates generated by Subroutine 3 Then we have
Amax (Anx(A)llss))) ll2
|
As
bll2 S
(k + 1)2
E.2 Implementation of SEP Oracle
We implement the SEP oracle in Definition 2 based on the classical Lanczos method with a random start, where the initial vector is chosen randomly and uniformly from the unit sphere (see, c.g.. 55 58). For completeness, the full algorithm is described in Subroutine [
To prove the correctness of our algorithm, first recall classical result in 51 the convergence behavior of the Lanczos method.
Proposition 3 (5I Theorem 4.2]). Consider symmetric matrix W and let A (W) and >a denote its largest and smallest eigenvalues, respectively. Then after k iterations of the Lanezos method with a random start, we find unit vectors u!!) and (dd) such that
P((W ( (W) <(> Xa(W))) I 1.648vde -R(2-)))) V7(26-1) P((Wu (d) (d)) Z (a(W) + <(>i(W) Aa(W))) - 1.648 (1)
As a corollary, to ensure that, with probability least q.
W)) and (Wu(d) (d) u (d)) An (W) (W) the number of iterations can be bounded by e 1/2 (11gg//1//g/// 7. Lemma 20. Let and S be the output of SEP(W; 6, g) in Subroutine
(1) ())>>>>>>)> ((W)
Then with probability at
least 1 9. they satisfy one of the following properties:
Case I: y S 1, then we have || [W Wlop < 1;
.
Case II: 7 > 1, then we have ////llop < 1, ||S]| 3 and (S, W B) 2 foror any
B such that |B B|opp S 1.
33

---

Subroutine 4 ;
Input: W s", > q (0, 1)
2: Set the number of iterations N1 + min log 11d .
3: Run Lanczos method with a random start for Ni iterations to get ()) and (d)) (cf. Proposition
Set Å < (Wul"), ((1) and da < (Wuld) (
5: Set Amax max -
6: if Amax S 1/2 then # Case 1:
which implies [ Wooop S
7: Return 7 = 2dmax and S =0 0
else Amax 2 then Case 11: y > and S defines a separating hyperplane A Aa then
9:
Return 7 2dmax and 3u (1) (1
11: else
12:
Return 2dmax and S = -3n (T)
13:
end if
14:
else
Amax <2
15:
Set the number of iterations N₂ min 4v28
log 11d +
16:
Run Lanczos method with random start for N2 iterations to get ( and ((d) (cf. Proposition 3
17:
Set λ < (Wa (()) and da <- (Wü(4) (()))
Set Amax max(Ai.- -da]
19:
if Amax 1 then
20:
Return 7 Amax 8 and S
21:
else if λι -Ad then
22:
Return Amas 8 and S ū (
23: else
24:
Return Amax + 8 and S ((d) ((()T
25:
end if
26:
end if
Proof. Note that in Subroutine 4 we first run the Lanczos method for e 1/2 log ltd
q² +
iterations,
where e 1. Thus, by Proposition with probability at least - q/2 we have
\ ( (2) (2)) z X1 (W) 1 (A)( (W))
(83)
& A (W) u(d) - (w(WW)++ + (AI(W)
Aa())).
(84)
Combining 83 and
84 we get
5 (Ai (w)) s - + Aw - Aa (W) s 2( - -d).
By plugging the above inequality back into 83 and 84 we further have
1
A1 (W) d + (A)(W) Aa (W)
Å + (d1
2d):
(85)
2
Aa (W) 2, (A (w)) > &x
(²1 Ed)-
1
(86)
2
Let Ämax max (1, Saj. By (85) and 86). can further bound the eigenvalues of W by
1 25max 2imax and >a(W) 2 -Smas 22muxx -22
(WW < Ámax +
2
Hence, that ||llop max(A (a())) 2imax- Now distinguish three cases.
(a) If Amax S 21
1 then we are in Case and the ExtEvec oracle outputs 7 = 22 S and
S = O. In this case, we indeed have W|] S 7 S 1.
34

---

(b) If Amax 2 2, then we are in Case II. In addition, if λy 2 - A4. then the ExtEvec oracle returns 2mmax and S = 3u [u (T) Similarly, - Å1. then the ExtEvec oracle returns y 2X and S = -3u () (u )T. Without loss of generality, consider the case where dι 2 -Às. Since || - 2 maxx 7, we have || /W/yllop S Also, since in is a unit vector, we have |||| = 3//u 11 1|2 3. Finally, for any B such that |Bllopp S 1, we have
(S, W B) = 3( Win
Bu Tuu 35max 3 2dmax = 7--1 1,
where we used the fact that Amax 2 the last inequality.
(c) If < Amax 2, we continue to run the Lanczos method for a total number of z - 1/2 log
ltd + 1 iterations, where fo. Thus, by Proposition 3 with probability at
q²t
least - q/2 we have
X1 A (Wa a()) 2 (W) 1 (()) - Aa(W)),
(87)
8
Aa A (Wi a S Ao (W) + 5 ((1(WW) - Aa(W)).
(88)
Let Amax max(A1, -ay- Since have \1( S 2dmax I 4 and Aa \w(W) > -2imax 2 -4, the above implies that X 2 (W) and da I (a( + 6. Hence, we can see that ||lpp = max(Ar ->a(W)) S Ammx + 6. We further consider two subcases.
(cl) If Amax < 1 - 6, then we are in Case I and the ExtEvec oracle outputs 7 Amax + and S = In this case, we indeed have || S S 1.
(c2) If Amax oracle returns y Amax
6, then we in Case II. In addition, x -Àa, then the ExtEvec
and S ())) Similarly, if -Xx Å₁ then the ExtEvec oracle returns 7 Amax + 6 and S (ū(dd) Without loss of generality, consider the case where Å 2 -Àa. Since || S Amax + 7. we have //W/ll I 1. Also, since (1) is a unit vector, we have |||| = 1/ =
Finally, for any B such that |[B|lop I 1, we have
(S, w - B) (6 T Wa(") (6() (o(D)TBB) Bi > 7max - 1 11 11---1111- -
This completes the proof.
E.3 Proof of Theorem 2
We divide the proof of Theorem into the following three lemmas.
Lemma 21. If we Algorithm 7 as specified in Theorem for N iterations, then the total number of line search steps can be bounded by 2N + logi/s (22). As a corollary, the total number of gradient queries is bounded by 3N. logi/a (abi)
02
Proof. In backtracking scheme, the number of steps in each iteration is given by logi/s (mp//t)))))))) Also note that 17 S fj// for all 2 Thus, we have
Uk +1 1ogi/8
N + logi/8 A0
M logi/8
f
0=0
E N + logi/9
+
n&k logi/9
+
+1)
110
M
k=
< 2N 1 logi/8
Furthermore, since îts 2 L2 for all k > 0, anrive at the conclusion.
35

---

Lemma 22.
The total number of matrix-vector product evaluations in the LinearSolver oracle is
bounded by Ne + Cuivool1 C12V
Li
where C11 and C12 absolute constants.
2
Proof. Our proof loosely follows the strategy in Subroutine before it terminates.
32
We first bound the number of steps required by
Lemma 23, Suppose A = I. Then Subroutine terminates after at most [Vatt Amax (A) - 1
iterations.
Proof.
Note that ||sl|2 2 ||s*||2 || 121 Also, since A > I, we have ||8k - -8"l12 S
|| - 5* s*)ll2 122 |xll2. Therefore, we have
1 - a||s
|rxll2 S al|s"122 -
|r I a |s 12
By using Lemma 19 we only need k
V
gttl Xmax (A) 1 to achieve ||Ask bl S a sk
Moreover, when the step size is smaller enough, we can show that Subroutine 3 will terminate in one iteration.
Lemma 24. Let A 1+ rBB. When n) S 2i, Algorithm
terminates in one iteration.
Proof. From the update rule of Subroutine 3
we can compute that Si
b Ab b, which implies TAb
|//22//2
||b]]
[b]
| || = |[b||
2
2
(A)/bb) A(A
Amax (A)
+ mL1
On the other hand, we also have
|rlll < ||b - bll Bb]| S n/Il/bll.
Moreover, when 7) S a 21, have 1/l1 < It.l. which implies that ||r || S a||si]]-
Now we upper bound the total number matrix-vecteon products in Algorithm 1 k-th iteration, the vector B&b at the beginning and reuse it to compute s when the step size 17+ < 2L, B'mbh1 01 it holds that
Note that at the
the LinearSolver oracle with A I+ B, where 7 7. We can store
aa And when
On +
1 1+ L+
S
3 MELL
a1
Thus, at the k-th iteration, the number matrix-vector products can be bounded by
MVA <1
Σ
a1 +
L1)
120m
1+ Σ
0n +
V 1'
320mm T2
01 + 2
S
MT..
01
Furthermore, we can bound that
N-1
N-1
N Σ
2(2 (3)
m S V
0a
Th
S V v0
vou
VB(1 - /B) ANN-
VB
k=
K=D
36

---

Note that {( [x*)=
No xpxx1 -00='112
21 Hence, have AN-1
1 xx'1
Thus, we can
2c
bound the total number of matrix-vector product evaluations by
MV = Σ
N -
MV & No
+ a1
01 2
(Voo1++ +
2(2 vB)
I Iz0
x=11
V3
v VB)
2c
k=0
L1 Zu - x21
Ne + CII 0OL1 + +12p
2E
where we define Ci1
+2 and C12 42
312-/3)
01
I-VB
I-VB
-3(1-
33-v33
Lemma 25. The total number of matrix-vector product evaluations in the SEP oracle is bounded by ON (log N.) log( vāN,
Proof. Note that we have N₂- 1v26;
log 1id in Subroutine where dN 1/(1 log(t= of
2)) and qi (2.5 1) log 2(t (t 1)). Thus, we have
T-1
T-1
+ 2) 0.25 log" 0. (t 2)
2.5v 44d(t + log (t+11
N = Σ
N₂
<
log
(89)
2v2
t=
t=0
VdN.
=
N. Vlog N. log
(90)
p
F Experiments
In our experiments, we consider the logistic regression problem and minimizing the log- Dog-sum-exp func- tion. Below we provide more details about the data generation scheme as well as the implementation of Nesterov's accelerated gradient method, BFGS, and our proposed A-QPNE algorithm.
Dataset generation. In the first experiment of logistic regression, the dataset consists of data points {(a, y: B)) where a, E R2 is the i-th feature vector and E {-1, 1} is its corresponding label. The labels {y: ) are generated by
Yi = sign( (ann(a))x)))
1,2, ...
where a. E Rd-1 and x" E Rd-1 are the underlying true feature vector and the underlying true parameter, respectively. Moreover, each entry of a and x" is drawn independently according the standard normal distribution (00,, 1). Note that the true feature vectors {a, |==1 are not given in our dataset; instead, we generate {a: by adding noises and appending an extra dimension {a{ 1==1-
Specifically, let a, |aţ n, 1] Rd where n, _N(0, I is the i.i.d. Gaussian noise vector and pd-1 denotes the all-one vector. In experiment, set 2, 000, d 150 and 0.8.
In the second experiment of log-sum-exp function, we follow a similar procedure as [6 to generate the dataset (a,,,b,)); b,)) where a, E IRd and b, R. First, we generate the auxiliary random vectors {ā. 1 by sampling each entry of á; uniformly and independently from the interval [-1, 1]. Moreover, we generate {b; independently from the standard normal distribution A(00, 1). Given {(ā,, b we define an auxiliary function f(x) = log(En l and finally let a; =
à v (0) for i = 1, ...., As discussed in [16), the purpose of this procedure is to ensure that V (0) = and thus is the unique minimizer f. In our experiment, we set n = d = 250.
NAG. We implemented a monotone variant of the Nesterov accelerated gradient method as described 5 Section 10.7.4]. Moreover, we determine the step size using a backtracking line search scheme.
BFGS. We implemented the classical BFGS algorithm, where the step size is determined by the line search scheme using implementation by Diane O' O'Leary!
http://www.cc umd. wttp://srr//oleerrrsssstsea
37

---

10
10
-O-NAG A--- -8 BFCS w/
LANEE
BEG3
109
10
[(ro) -(e')
1
660
10
00000
10"
10:10
0
00966
15"5
16
0
0.5
1.5
0
5
10
15
20
Time/a
Thooo/
(a) Logistic regression
(b) Log-sum-exp
Figure Plots of Suboptimality gap in terms of the running time.
A-QPNE (our method). We implemented proposed A-QPNE method following the pseudocode in Algorithm[] | where the line search scheme is given in Subroutine and the Hessian approximation update is given in Subroutine Moreover, the implementations of the LinearSolver oracle and the SEP oracle are given by Subroutines 3 and respectively.
F1 Additional Plots
In Fig. compare the performance of our proposed A-QNPE method with NAG and BFGS in terms of the runming time. All experiments are conducted using MATLAB R2021b on a MacBook Pro with Apple M1 chip and 16GB RAM. We observe from Fig. (3aa that method requires less running time than NAG due to its faster convergence, especially when we are seeking a solution of high accuracy. On the other hand, there cases where method is slower than NAG in terms of the running time, as shown in Fig. b). This is because, in this case, the cost gradient computation is comparable to the cost of matrix-vector product evaluation, and therefore method incurs higher computational cost per iteration than NAG.
38


Reviews:
[{"review_id": "t8GpLKUv4j", "reviewer": "Reviewer_MQHh", "summary": "This paper uses the optimal and adaptive Monteiro Svaiter acceleration framework to create a quasi-Newton method that solves unconstrainted convex problems with Lipschitz gradients and Lipschitz hessians at Nesterov's accelerated rate $O(1/k^2)$ but when the number of iterations is greater enough than the dimension, namely $\\Omega(d\\log d)$ it provides better convergence guarantees in terms of gradient oracle complexity. The number of operations is superlinear $\\widetilde{O}(N_\\epsilon^{1.25})$ in the number of gradient queries $N_\\epsilon$ and memory that is quadratic in the dimension is used.\n", "questions": ".", "limitations": "see above", "rating": 7, "confidence": 4, "soundness": 4, "presentation": 3, "contribution": 4, "strengths": "Authors use a great range of technical and powerful tools like optimal and adaptive MS acceleration, they use projection free online learning with a separation oracle, the conjugate gradients method and Lanczos algorithms.  The results are novel and interesting. \n", "weaknesses": "The abstract (and some parts of the paper) says that the method matches the optimal rate O(1/k^2) by NAG, and this rate is known optimal for functions that are convex smooth. Here your setting is convex smooth + Lipchitz Hessians. It would be good to explicitly comment on how that construction has a hessian Lipschitzness L_2 = 0 (since it is a quadratic) and therefore it also applies to your setting.\n\nI guess you can modify your algorithm to get results under the additional assumption of \\mu strong convexity. This could take a lot of work, but on the other hand you can follow the spirit of the usual reductions to show rates for your algorithm under strong convexity by using a sequence of restarts. So I would suggest to add a remark with this. The argument would be like this, you run the algorithm in stages and after each stage you guarantee you halve the distance to x^ast and so you only need to run a logarithmic number of stages. In order to do that, you want to run the algorithm for k iterations such that you guarantee the last equality here:\n$$\n\\mu/2 \\|x_t-x^\\ast\\|^2 \\leq f(x_t) - f(x^\\ast) \\leq O(\\text{your bound}) = \\mu/8 \\|x_0-x^\\ast\\|^2\n$$\nThe first inequality is strong convexity and the second one is your guarantee. In the regime in which \\mu is small enough, the number of iterations necessary in each stage should be \\Omega(d\\log(d)) and so your improved bound kicks in.\nAfter doing this, it would be desirable to have a comparison of these results with the quasi-Newton results mentioned in the introduction that apply to strongly convex functions.\n\n\nSimilarly, regarding \"However, all of the results above only apply under the restrictive assumption that the objective function f is strictly or strongly convex. In the more general setting where f is merely convex...\". Given the reduction that solves convex smooth  problems by regularizing with + \\epsilon/ R^2 \\norm{x_0-x^\\ast}^2, where R is an upper bound on the initial distance to the minimizer, you should discuss what one can get with previous methods under the reduction. That is, the strongly convex setting is not necessarily restrictive and a comparison / discussion should be done.\n\nL839 \"subourtine\" -> \"subroutine\" \n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the insightful comments. We address your concerns below.\n\n---\n\n**Q1 In your setting the function is convex smooth + Lipchitz Hessians. It would be good to comment on how the worst-case construction also applies to your setting.** \n\n**A1** This is an excellent point. The lower bound of $\\Omega(1/k^2)$ for the class of convex smooth functions is established by a worst-case quadratic function, whose Hessian is constant (Lipschitz continuous with $L_2 = 0$). Therefore, the additional assumption of Lipschitz Hessian does not eliminate this worst-case construction from the considered problem class, and thus the $\\Omega(1/k^2)$ lower bound also applies to our setting. Thanks for raising this excellent point and we will add it to the revision. \n\n---\n\n**Q2 Use restart to extend the results to the strongly convex setting.** \n\n**A2** This is a very good suggestion. Indeed, it is possible to use the restarting technique to further extend our result to the strongly convex setting, as suggested by the reviewer. However, it appears that such a reduction would only lead to a linear rate with a better dependence on the condition number, instead of a superlinear rate that we would expect from a quasi-Newton method in the strongly convex setting.\n\nTo be more precise, if we follow the arguments of the reviewer and run our algorithm in multiple stages, where we guarantee that the distance to $x^\\*$ is halved after each stage, we obtain from Theorem 1 that \n$$f(x_t) - f(x^\\*) = O\\left( \\min\\left\\\\{\\frac{L_1\\\\|x_0-x^\\*\\\\|^2}{t^2}, \\frac{L_1\\sqrt{d}\\\\|x_0-x^\\*\\\\|^2}{t^{2.5}}\\right\\\\} \\right).$$ To get this bound, we upper bound $\\\\|B_0 - \\nabla^2 f(z_0)\\\\|_F^2$ by $L_1^2 d$, ignore the log factor in (19) and only focus on the dominant term for simplicity. \nAs the reviewer points out, by using strong convexity we have $\\\\|x_t-x^\\*\\\\| \\leq \\frac{1}{2}\\\\|x_0-x\\^\\*\\\\|$ if $f(x_t)-f(x^\\*) \\leq \\frac{\\mu}{8}\\\\|x_0-x^\\*\\\\|^2$. Hence, the number of iterations required in each stage can be bounded by $O(\\min\\\\{({\\frac{L_1}{\\mu}})^{0.5}, d^{0.2}(\\frac{L_1}{\\mu})^{0.4}\\\\})$, which implies a total complexity of $$ O\\left(\\min\\left\\\\{\\left(\\frac{L_1}{\\mu}\\right)^{0.5}, d^{0.2}\\left(\\frac{L_1}{\\mu}\\right)^{0.4}\\right\\\\} \\log \\frac{1}{\\epsilon}\\right).$$ \n\nIn the regime where $d \\leq \\sqrt{\\frac{L_1}{\\mu}}$, the obtained complexity bound outperforms NAG in terms of the dependence on the condition number. \nOn the other hand, we note that several papers [1,2] have established a local non-asymptotic superlinear rate of the form $O((1/\\sqrt{k})^k)$ for classical quasi-Newton methods and their variants. More recently, a global non-asymptotic superlinear rate is also shown for a quasi-Newton proximal extragradient method [3]. In comparison, the restarting scheme described above can only achieve global linear convergence.  \n\nWhile the main focus of this paper is on demonstrating a provable gain for a quasi-Newton-type\nmethod over NAG in the convex setting, the above argument for extending our algorithm to the strongly convex setting is an interesting observation and we will add it as a remark to our revised paper. \n\n---\n\n**Q3 Discuss what one can get with previous methods under the regularization reduction.**  \n\n**A3** \nThanks for your insightful comment. As the reviewer rightly pointed out, one can regularize $f$ with $\\frac{\\epsilon}{R^2}\\\\|x-x_0\\\\|^2$ to reduce a convex smooth problem into a strongly-convex one with $\\mu = \\frac{\\epsilon}{R^2}$. However, to the best of our knowledge, applying this reduction directly to the existing analysis of quasi-Newton methods would not lead to a global complexity bound better than the one for NAG, as we elaborate next.\n\n- The results in [1,2] are crucially based on local analysis and require the initial point $x_0$ to be close enough to the optimal solution $x^*$. However, it is unclear how to explicitly bound the number of iterations before the iterate enters the local neighborhood, and even if this can be done, it seems unlikely that the total complexity would be better than NAG, as these results only provide a local convergence analysis. \n\n- The result from [3] seems to be the only strongly-convex result that can be compared with NAG in the convex setting using the mentioned regularization idea, as it provides a global convergence analysis with an explicit overall complexity bound. Based on the discussions in Appendix D.2 of [3], the authors showed a global complexity bound in the form of \n$O\\left(\\min\\left\\\\{\\frac{L_1}{\\mu} \\log \\frac{1}{\\epsilon},d^{\\frac{1}{3}}\\left(\\frac{L_1}{\\mu}\\log\\frac{1}{\\epsilon}\\right)^{\\frac{2}{3}} \\right\\\\}\\right)$ for strongly-convex objectives. Since we have $\\mu = \\frac{\\epsilon}{2D^2}$ under the reduction, this translates into a complexity bound of $\\tilde{O}\\left(\\min\\left\\\\{ \\frac{1}{\\epsilon}, d^{\\frac{1}{3}}\\left(\\frac{1}{\\epsilon}\\right)^{\\frac{2}{3}}\\right\\\\}\\right)$ for convex problems, which is worse than the bound $O(\\left(\\frac{1}{\\epsilon}\\right)^{\\frac{1}{2}})$ by NAG. This is conceivable since there is no form of acceleration in the proposed method of [3].  \n\nThus, we conclude that simply applying the standard reduction to the existing analysis would not result in a complexity bound better than NAG. This highlights the need for a distinct algorithm and analysis specifically tailored to the convex setting, as presented in this paper. We will add a remark to the revised paper regarding this point.\n\n--- \n\n**Typo.**\nThanks for catching the typo. We will fix this in the revision. \n\n---\n\nReferences:\n\n[1] Q. Jin and A. Mokhtari. Non-asymptotic superlinear convergence of standard quasi-Newton methods, 2022.\n\n[2] A. Rodomanov and Y. Nesterov. New results on superlinear convergence of classical quasi-Newton methods, 2021.\n\n[3] R. Jiang, Q. Jin, and A. Mokhtari. Online learning guided curvature approximation: A quasi-Newton method with global non-asymptotic superlinear convergence, 2023."}, {"review_id": "XXGt5CVsT5", "reviewer": "Reviewer_CWLw", "summary": "This paper proposes an accelerated quasi-Newton proximal extragradient method for solving unconstrained smooth convex optimization problems. The algorithm  can achieve a convergence rate of $\\mathcal{O}\\bigl(\\min\\{\\frac{1}{k^2}, \\frac{\\sqrt{d\\log k}}{k^{2.5}}\\}\\bigr)$, where $d$ is the problem dimension and $k$ is the number of iterations. In particular, in the regime where $k = \\mathcal{O}(d)$, our method matches the optimal rate of $\\mathcal{O}(\\frac{1}{k^2})$ by Nesterov's accelerated gradient (NAG). Moreover, in the the regime where $k = \\Omega(d \\log d)$, it outperforms NAG and converges at a faster rate of $\\mathcal{O}\\bigl(\\frac{\\sqrt{d\\log k}}{k^{2.5}}\\bigr)$. ", "questions": "No", "limitations": "No", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "This paper proposes an accelerated quasi-Newton proximal extragradient method for solving unconstrained smooth convex optimization problems. The algorithm  can achieve a convergence rate of $\\mathcal{O}\\bigl(\\min\\{\\frac{1}{k^2}, \\frac{\\sqrt{d\\log k}}{k^{2.5}}\\}\\bigr)$, where $d$ is the problem dimension and $k$ is the number of iterations. In particular, in the regime where $k = \\mathcal{O}(d)$, our method matches the optimal rate of $\\mathcal{O}(\\frac{1}{k^2})$ by Nesterov's accelerated gradient (NAG). Moreover, in the the regime where $k = \\Omega(d \\log d)$, it outperforms NAG and converges at a faster rate of $\\mathcal{O}\\bigl(\\frac{\\sqrt{d\\log k}}{k^{2.5}}\\bigr)$. ", "weaknesses": "I have a concern about the total computation cost. \nLet us consider the finite-sum form, that is, $f(x) = \\frac{1}{n} \\sum_{i=1}^nf_i(x)$ with $n = \\mathcal{O}(d)$.\nThen, the algorithm takes $N_\\eps$ gradient queries which implies $\\mathcal{O}(N_\\eps d^2)$ computation cost.\nBy Theorem 2 (c), the cost of computing SEP oracles, the computation cost is $ \\mathcal{O}(N_\\eps^{1.25} d^2) $. \nIn this case, A-QNPE does not have advantages over  Nesterov's accelerated gradient (NAG).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their comment. We address your concern below. \n\n---\n\n**Q1 Comparison with NAG in terms of the total computational cost.** \n\n**A1** Thanks for raising this point. For easier comparison, we follow the suggestion by Reviewer H1ti and summarize the computation cost of NAG and our proposed method A-QNPE to achieve an $\\epsilon$ accuracy in the following table.  \n\n|        | Gradient queries  | Matrix-vector products\n------- |    -------------  | -------------\nNAG     | $O(\\frac{1}{\\epsilon^{0.5}})$ | N.A.\nA-QNPE  | $\\tilde{O}(\\min\\\\{\\frac{1}{{\\epsilon}^{0.5}},\\frac{d^{0.2}}{\\epsilon^{0.4}}\\\\})$  | $\\tilde{O}(\\min\\\\{\\frac{d^{0.25}}{\\epsilon^{0.5}}, \\frac{1}{\\epsilon^{0.625}}\\\\})$\n\n\nWe observe that A-QNPE outperforms NAG in terms of gradient query complexity: it makes fewer or equal gradient queries, especially when $\\epsilon < \\frac{1}{d^2}$. On the other hand, A-QNPE requires additional matrix-vector product computations to implement the LinearSolver and SEP oracles. \n\n\nTo assess the overall computation cost, we have to consider the cost of gradient computation, which varies depending on the specific problems. As a concrete example, \nconsider the finite-sum minimization problem $f(x) = \\frac{1}{n} \\sum_{i=1}^nf_i(x)$. In this case, one gradient query typically costs $O(nd)$, while one matrix-vector product costs $O(d^2)$. Thus, the total computation cost of NAG and A-QNPE can be bounded by $O(\\frac{nd}{\\epsilon^{0.5}})$ and $ O(\\frac{n d^{1.2}}{\\epsilon^{0.4}}+ \\frac{d^{2.25}}{\\epsilon^{0.5}})$, respectively. In particular, our method will incur a lower computation cost when $\\epsilon \\ll \\frac{1}{d^2}$ and $n \\gg d^{1.25}$.\n\nAs a final remark, we acknowledge that our method may be faster or slower than NAG, depending on the specific problem. Nevertheless, we would like to highlight that this is the first work to **theoretically demonstrate that quasi-Newton-type methods can outperform NAG in certain regimes**. Indeed, previous works [1,2] on quasi-Newton methods provide a convergence rate matching NAG, and as a result their overall computational cost will always be larger than NAG in theory. Thus, we believe our paper is an important conceptual advance for quasi-Newton methods and we leave the task of further reducing the computation cost as future work.     \n\n---\n\nReferences: \n\n[1] K. Scheinberg and X. Tang. Practical inexact proximal quasi-Newton method with global complexity analysis, 2016\n\n[2] H. Ghanbari and K. Scheinberg. Proximal quasi-Newton methods for regularized convex optimization with linear and accelerated sublinear convergence rates, 2018."}, {"review_id": "gZGhNCdD41", "reviewer": "Reviewer_H1ti", "summary": "This paper proposed a novel quasi-Newton method with faster global convergence rate. The algorithm uses the framework of MS acceleration and updates the Hessian estimator via online learning. The obtained convergence rate is impressive, it firstly show a faster global rate of quasi-Newton method which cannot be achieved by the first order methods. ", "questions": "1. Can you compare the detailed computation cost in a table to make the results more clear?\n2. Can you conduct some experimental results in terms of the running time?\n", "limitations": "N/A", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The proposed algorithm is novel. The global convergence rate of $\\tilde{\\mathcal{O}}(\\sqrt{d}/k^{-2.5})$ is a significant theoretical result. ", "weaknesses": "1. Although the iteration complexity is impressive and significant, the total computation complexity of the proposed method is not very satisfactory. Compared with the first-order method, it requires additional Hessian-vector products with the complexity of $\\min\\lbrace d^{0.25}/\\epsilon^{0.5},1/\\epsilon^{0.625}\\rbrace$ which may lead to an even higher computation cost than NAG.\n\n2. The proposed method is complicated and may not practical to use. For the experimental part, the authors do not present the comparison between their method and the baselines in term of running time.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their comment. We address your concerns below. \n\n---\n\n**Q1 Compare the detailed computation cost in a table.** \n\n**A1** \nThanks for your suggestion. In the following table, we summarize the detailed computation cost of NAG and our method A-QNPE, and we will also include it in our revision. In addition, we would like to clarify that our method requires extra matrix-vector products, rather than Hessian-vector products. Indeed, as a quasi-Newton method, we only need access to the gradient oracle.  \n\n|        | Gradient queries  | Matrix-vector products\n------- |    -------------  | -------------\nNAG     | $O(\\frac{1}{\\epsilon^{0.5}})$ | N.A.\nA-QNPE  | $\\tilde{O}(\\min\\\\{\\frac{1}{{\\epsilon}^{0.5}},\\frac{d^{0.2}}{\\epsilon^{0.4}}\\\\})$  | $\\tilde{O}(\\min\\\\{\\frac{d^{0.25}}{\\epsilon^{0.5}}, \\frac{1}{\\epsilon^{0.625}}\\\\})$\n\n\n\n\nWe observe that A-QNPE outperforms NAG in terms of gradient query complexity: it makes fewer or equal gradient queries, especially when $\\epsilon < \\frac{1}{d^2}$. On the other hand, A-QNPE requires additional matrix-vector product computations to implement the LinearSolver and SEP oracles. \n\n\nTo assess the overall computation cost, we have to consider the cost of gradient computation, which varies depending on the specific problems. As a concrete example, \nconsider the finite-sum minimization problem $f(x) = \\frac{1}{n} \\sum_{i=1}^nf_i(x)$. In this case, one gradient query typically costs $O(nd)$, while one matrix-vector product costs $O(d^2)$. Thus, the total computation cost of NAG and A-QNPE can be bounded by $O(\\frac{nd}{\\epsilon^{0.5}})$ and $ O(\\frac{n d^{1.2}}{\\epsilon^{0.4}}+ \\frac{d^{2.25}}{\\epsilon^{0.5}})$, respectively. In particular, our method will incur a lower computation cost when $\\epsilon \\ll \\frac{1}{d^2}$ and $n \\gg d^{1.25}$.\n\nAs a final remark, we acknowledge that our method may be faster or slower than NAG, depending on the specific problem. Nevertheless, we would like to highlight that this is the first work to **theoretically demonstrate that quasi-Newton-type methods can outperform NAG in certain regimes**. Indeed, previous works [1,2] on quasi-Newton methods provide a convergence rate matching NAG, and as a result their overall computational cost will always be larger than NAG in theory. Thus, we believe our paper is an important conceptual advance for quasi-Newton methods and we leave the task of further reducing the computation cost as future work. \n\n---\n\n**Q2 Conduct some experimental results in terms of the running time?**\n\n**A2**\nThanks for your suggestion. We have included additional plots in terms of the running time; please check Figs. 1 and 2 in the attached pdf. All experiments are conducted using MATLAB R2021b on a MacBook Pro with an Apple M1 chip and 16GB RAM. \n\nSpecifically, in Fig. 1 of the attached pdf file, we consider the logistic regression problem $f(x)= \\frac{1}{n}\\sum_{j=1}^n \\log(1+e^{-y_j \\langle a_j, x\\rangle})$ as described in the paper, where the dimension is $d = 150$ and the number of samples is $n=2000$. As shown in Fig. 1(c), if we are seeking a solution of high accuracy,  our method can require less running time than NAG due to its faster convergence.\n\nIn addition, in Fig. 2 of the attached pdf file, we consider the log-sum-exp function following the suggestion of Reviewer N8Fk. The loss function is given by $f(x) = \\log( \\sum_{j=1}^n e^{\\langle a_j, x\\rangle - b_j}) $, where the dimension is $d=150$, the number of samples is $n=150$, and we follow the procedure in [1] to generate $\\\\{a_j\\\\}$ and $\\\\{b_j\\\\}$. In this case, Fig. 2(c) shows that the run-time performance of A-QNPE is comparable to that of NAG. \n\nFinally, we would like to emphasize that the primary objective of our numerical experiments is to validate our theoretical discovery that A-QNPE can attain a faster convergence rate than NAG. With a more meticulous implementation, there is potential to enhance the practical efficacy of our method, and we have deferred this as a prospect for future investigation.\n\n\n[1] A. Rodomanov and Y. Nesterov. Greedy quasi-Newton methods with explicit superlinear convergence, 2021"}, {"review_id": "hG6IGmr9AE", "reviewer": "Reviewer_N8Fk", "summary": "This paper proposes a first-order optimization method for convex optimization that requires gradient computations and matrix-vector products. Based on recent advances on quasi-newton methods, this method converges at a rate that matches the optimal rate ($1/k^2$) when $k = \\Omega(d)$ and improves upon the optimal rate when $k >>d$. In order to achieve this rate the authors consider a quasi-Newton method with a backtracking line search and a projection-free online algorithm to approximate the hessian of the objective. ", "questions": "Can you develop why you cannot use standard projection-free no-regret methods? such as\n\n\nNew Projection-free Algorithms for Online Convex Optimization with Adaptive Regret Guarantees \nDan Garber and Ben Kretzu\n\n\nI assume it would not be straightforward since your notion of regret is different (dynamic regret) but I would be curious to know how your method compares with respect to this related work.\n", "limitations": "They authors mostly do address the limitation of their work. \n\nThey mention that their algorithm requires many matrix-vector products and, thus, their rate is faster than NAG if the gradient computation is considered the bottleneck. \n\nThey mention that their algorithm is slower than BFGS in practice. ", "rating": 8, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "This paper has several strengths:\n- the rate proposed improved upon the best-known rates. \nIn particular, it is faster than the convergence rate of NAG when the number of iterations is larger than the dimension. \n- The authors really try to give the main steps of the results in the main paper. \n- While the paper is quite technical, I found that the contributions are well-presented, \n", "weaknesses": "I found the paper could be slightly improved:\n- The experiments could support the theory more:\n  - The logistic loss is strictly convex. Thus it seems that the rate in the experiments should be super-linear. It is not clear in the experiments if we are observing these super-linear rates or if we see sublinear convergence. I would instead try the log-sum exp loss (logistic regression for multi-class), which is not strictly convex \n  - I would be curious to see if, actually, one can find a practical loss for which one can observe the sublinear rate $1/k^{2.5}$.\n  - A plot with time as a x-axis could illustrate the fact that in many situations, matrix-vector product computations are not the bottleneck.\n- Some details on the Online method are missing in the main paper:\n   - The subroutine could be described (at least intuitively) \n   - The main result about the regret bound could be stated (basically logarithmic regret) \n   - In theorem 1 the probability p does not appear in the bound. I understand it only appears in the matrix-vector product complexity but it makes the statement quite odd. I would suggest the authors to at least mention it in a remark. ", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the detailed comments. We address your concerns below. \n\n---\n\n**Q1 It is not clear in the experiments if we are observing super-linear rates or if we see sublinear convergence. Try the log-sum exp loss, which is not strictly convex.** \n\n**A1** \nThis is a very good observation. We conducted new experiments on the log-sum-exp function $f(x) = \\log( \\sum_{j=1}^n e^{\\langle a_j, x\\rangle - b_j}) $, where the dimension is $d=150$, the number of samples is $n=150$, and we follow the procedure in [1] to generate $\\\\{a_j\\\\}$ and $\\\\{b_j\\\\}$. As shown in Fig. 3 in the attached pdf, we observe that both NAG and A-QNPE converge at a sublinear rate. \n\n---\n\n**Q2 I would be curious to see if one can find a practical loss for which one can observe the sublinear rate $O(1/k^{2.5})$.**\n\n**A2**\nThanks for raising this point. So far we have not observed the sublinear rate $O(1/k^{2.5})$ in all of our experiments. Indeed, our theoretical result is pessimistic by nature, and the empirical performance of A-QNPE can be better than what the theory predicts. For instance, in our new experiment, Fig. 3(c) indicates that on this specific problem NAG converges at the rate of $O(1/k^{3})$, while A-QNPE converges at the rate of $O(1/k^5)$. \n\n--- \n\n**Q3 A plot with time as a x-axis.** \n\n**A3** \nThanks for your suggestion. We have included additional plots in terms of the running time in the attached pdf.  All experiments are conducted using MATLAB R2021b on a MacBook Pro with an Apple M1 chip and 16GB RAM. \n\nSpecifically, in Fig. 1 of the attached pdf file, we consider the logistic regression problem $f(x)= \\frac{1}{n}\\sum_{j=1}^n \\log(1+e^{-y_j \\langle a_j, x\\rangle})$ as described in the paper, where the dimension is $d = 150$ and the number of samples is $n=2000$. As shown in Fig. 1(c), if we are seeking a solution of high accuracy,  our method can require less running time than NAG due to its faster convergence.\n\nIn addition, in Fig. 2 of the attached pdf file, we consider the log-sum-exp function $f(x) = \\log( \\sum_{j=1}^n e^{\\langle a_j, x\\rangle - b_j}) $, where $d=150$, $n=150$, and we follow the procedure in [1] to generate $\\\\{a_j\\\\}$ and $\\\\{b_j\\\\}$. In this case, Fig. 2(c) shows that the run-time performance of A-QNPE is comparable to that of NAG. \n\nFinally, we would like to emphasize that the primary objective of our numerical experiments is to validate our theoretical discovery that A-QNPE can attain a faster convergence rate than NAG. With a more meticulous implementation, there is potential to enhance the practical efficacy of our method, and we have deferred this for future investigation.\n\n--- \n\n**Q4 Some details on the Online method are missing in the main paper.** \n\n**A4**\nThanks for raising this point. Due to the page limit, we have to relegate the details of our online learning algorithm to the appendix (Sections C.3 and C.4). In the revision, we will describe the subroutine in more detail and present the key regret bound (Lemma 12 in Section D) in the main paper. \n\n---\n\n**Q5 In theorem 1 the probability $p$ does not appear in the bound.**\n\n**A5** Thanks for raising this point. Note that our algorithm relies on the SEP oracle in Definition 2, which succeeds with a certain probability. To prove Theorem 1, we first show that, with probability at least $1-p$, every call of the SEP oracle is successful during the execution of our algorithm. Conditioned on this event, we proceed to prove our convergence rates, and thus they do not depend on $p$.  Following your suggestion, we will add a remark on this in our revision. \n\n---\n\n**Q6 Can you develop why you cannot use standard projection-free no-regret methods?** \n\n**A6**\nThis is a good point. We note that standard projection-free no-regret methods, such as online Frank-Wolfe [2], are based on a linear minimization oracle (LMO). Unfortunately, implementing the LMO in our setting is also computationally expensive. Specifically, the constraint set in our online learning problem is given by $\\mathcal{Z} = \\\\{\\mathbf{B} \\in \\mathbb{S}^d_+: 0 \\preceq \\mathbf{B} \\preceq L_1 \\mathbf{I}\\\\}$. Consider the linear minimization problem $\\min_{\\mathbf{X} \\in \\mathcal{Z}} \\langle \\mathbf{A}, \\mathbf{X}  \\rangle.$\nWe first need to compute the eigendecomposition $\\mathbf{A} = \\mathbf{V}\\mathbf{\\Lambda} \\mathbf{V}^\\top$, where $\\mathbf{V}$ is an orthogonal matrix and $\\mathbf{\\Lambda} = \\mathrm{diag}(\\lambda_1,\\dots,\\lambda_d)$ is a diagonal matrix. Then the solution is given by $\\mathbf{V}\\mathbf{\\Lambda}' \\mathbf{V}^\\top$, where  $$\\lambda_k' = \\begin{cases}\n    0, & \\text{if } \\lambda_k \\geq 0; \\\\\\\\\n    L_1, & \\text{otherwise.}\n\\end{cases}$$\nHence, implementing the LMO  would require a full matrix eigendecomposition, which requires $O(d^3)$ arithmetic operations in general.  \n\nOn the other hand, as we discuss in this paper, the (approximate) separation oracle of $\\mathcal{Z}$ is more efficient to implement, since it only requires computing the two extreme eigenvectors and eigenvalues of the given matrix. To the best of our knowledge, only the two recent papers [3,4] consider no-regret algorithms based on the separation oracle. In this paper, we develop our online learning algorithm based on [4] for two main reasons: (1) the algorithm in [4] appears to be simpler as it only requires one call to the separation oracle per iteration; (2) it is relatively straightforward to allow inexactness of the separation oracle. That said, we think it might also be possible to adapt the algorithm in [3] and achieve similar convergence guarantees. \n \n\n---\n\n[1] A. Rodomanov and Y. Nesterov. Greedy quasi-Newton methods with explicit superlinear convergence, 2021.\n\n[2] E. Hazan and S. Kale. Projection-free online learning, 2012.\n\n[3] D. Garber and B. Kretzu. New Projection-free Algorithms for Online Convex Optimization with Adaptive Regret Guarantees, 2022. \n\n[4] Z. Mhammedi. Efficient projection-free online convex optimization with membership oracle, 2022."}]

Author rebuttal:
We thank all reviewers for their time and effort in evaluating our paper. Following the suggestions by **Reviewer H1ti** and **Reviewer N8Fk**, we have included additional plots in the attached pdf file. 

- In Fig. 1, we consider the logistic regression problem $f(x)= \frac{1}{n}\sum_{j=1}^n \log(1+e^{-y_j \langle a_j, x\rangle})$ as described in the paper, where the dimension is $d = 150$ and the number of samples is $n=2000$. As shown in Fig. 1(c), if we are seeking a solution of high accuracy,  our method can require less running time than NAG due to its faster convergence.

- In Fig. 2, we consider the log-sum-exp function following the suggestion of **Reviewer N8Fk**. The loss function is given by $f(x) = \log( \sum_{j=1}^n e^{\langle a_j, x\rangle - b_j}) $, where the dimension is $d=150$, the number of samples is $n=150$, and we follow the procedure in [1] to generate $\\{a_j\\}$ and $\\{b_j\\}$. In this case, Fig. 2(c) shows that the run-time performance of A-QNPE is comparable to that of NAG. 

- In Fig. 3, we plot both the suboptimality gap $f(x_k) - f(x^*)$ and the number of iterations on a log scale for the log-sum-exp experiment. For this specific problem, we can observe empirically that NAG converges at a sublinear rate of $O(1/k^{3})$, while A-QNPE converges at a faster rate of $O(1/k^5)$. 

---

[1] A. Rodomanov and Y. Nesterov. Greedy quasi-Newton methods with explicit superlinear convergence, 2021

Meta-review: The reviewers reached a consensus about accepting this paper, which two of the reviewers expressing excitement about the progress this paper makes in theoretically establishing quasi-Newton methods. Consequently, I recommend acceptance of the paper.

When making the camera-ready revision, please take care to implement all the changes proposed during the review and rebuttal process, and particularly the ones concerning the discussion and empirical evaluation of the computational cost of the proposed method.
Decision: Accept (spotlight)