Blake Woodworth

Konstantin Mishchenko

Francis Bach

Two Losses Are Better Than One: Faster Optimization Using a Cheaper Proxy

We present an algorithm for minimizing an objective with hard-to-compute gradients by using a related, easier-to-access function as a proxy. Our algorithm is based on approximate proximalpoint iterations on the proxy combined with relatively few stochastic gradients from the objective. When the difference between the objective and the proxy is δ-smooth, our algorithm guarantees convergence at a rate matching stochastic gradient descent on a δ-smooth objective, which can lead to substantially better sample efficiency. Our algorithm has many potential applications in machine learning, and provides a principled means of leveraging synthetic data, physics simulators, mixed public and private data, and more.

Introduction

In many machine learning problems, it is difficult to access the objective function, e.g., to compute its gradients for use in training. This difficulty can come from a wide range of sources. For example, it might be costly or time consuming to collect labelled training examples, leading to a small training set. In some applications, computing the gradient of the loss might require waiting for an actual robot to execute a policy in a physical environment, which could be very slow. Moreover, physical execution can risk damage to the robot or its surroundings. Finally, individuals' data may be subject to privacy or legal constraints, limiting access for gradient computations. In all of these cases, we face challenges that complicate training, and can lead to poor performance using off-the-shelf algorithms like stochastic gradient descent (SGD).

A natural solution to the problems above is to find a second, easier-to-access "proxy" function that can be used in place of the real objective. When few samples are available, syn-thetic data could be used to approximate the objective, or in applications with sensitive data, we could train a generative model using a privacy-preserving algorithm and define the proxy using synthetic samples. To mitigate the cost and risk of operating a physical robot, one could use a simulator that is faster, cheaper, and safer.

Unfortunately, directly optimizing a hand-crafted proxy objective might not improve performance on the original problem, because the surrogate loss can introduce a bias. We aim to remedy this in our work, and our main contribution is a general algorithm that exploits curvature information from the proxy loss in order to better optimize the objective that we care about.

The problem

Our ultimate aim is to solve unconstrained convex optimization problems of the form min w∈R d L(w) .

(1)

Here, L could represent the population loss of a machine learning model parametrized by weights w; it could also be the training loss, or any other objective function to be minimized. To solve (1), we use a proxy function L and stochastic gradients g k approximating ∇L(w k ) and update w k+1 ≈ arg min w L(w)

Proxy loss

+ g k -∇ L(w k ), w -w k Bias correction + 1 2η w -w k 2 Regularization . (2)
Thus, each iteration uses just a single stochastic gradient from L, but requires a heavier computation based on L. If we take L ≡ 0 , we note that the update ( 2) is equivalent to one stochastic gradient descent step

w k+1 = w k -ηg k . (3)
Therefore, SGD is a natural point of comparison for our algorithm, and the main question is to what extent updating using (2) (with non-zero L) rather than (3) is worth the additional computational cost. We will show that when L is sufficiently "similar" to L, our algorithm's updates, (2), can converge to the minimum of L in substantially fewer iterations than would be needed for SGD, (3). Accordingly, the number of stochastic gradients from L that our algorithm needs is less than what would be needed by SGD.

For this to work, we obviously require that L somehow resembles L, and the particular notion of similarity that we consider is that the function h := L -L is differentiable and has δ-Lipschitz continuous gradients, which, in the case that L and L are twice-differentiable, is equivalent to requiring that ∀ w ∇ 2 L(w) -∇ 2 L(w) op ≤ δ. This measure, sometimes referred to as "δ-Hessian similarity", is common in the optimization literature (see, e.g., Mairal, 2013a;2015;[START_REF] Arjevani | Communication complexity of distributed convex learning and optimization[END_REF][START_REF] Kovalev | Optimal gradient sliding and its application to distributed optimization under similarity[END_REF][START_REF] Chayti | Optimization with access to auxiliary information[END_REF]. For certain problems, it can be easy to construct a proxy that satisfies this similarity condition for small δ. For instance, the Hessian of a leastsquares objective is simply the feature covariance matrix, so an appropriate proxy can be defined using any estimate of the covariance matrix that is δ-accurate in the operator norm. In Section 3.1 we further show that for least squares and logistic regression, proxies with δ = 0 can be constructed without using any labels.

Our contributions

We present Algorithm 1, which substitutes cheaper accesses to a proxy loss in place of costly accesses to the objective itself. In Section 2.2, we prove convergence guarantees for our algorithm in the convex and strongly convex settings, and in Section 2.4, we show that these guarantees imply statistical optimality as well as a dependence on the similarity parameter, δ, in the place of other properties of L itself, such as its potentially much larger smoothness parameter. On a technical level, our analysis improves over many similar methods by allowing each iteration's proximal-point subproblem to be solved inexactly, with the inexactness captured by a simple criterion that can be evaluated at the time of execution. Furthermore, in Section 2.3 we show that this criterion can be satisfied efficiently owing to the strong convexity of the proximal-point subproblem's objective. In Section 3, we discuss potential applications of our algorithm where the target objective is costly to access and a suitable proxy is available. In Section 4, we conduct experiments to show the efficacy of our algorithm in realistic problems, including one with a non-convex objective. Finally, in Section 5, we discuss several possible extensions of our work and we provide a preliminary analysis of our algorithm for the non-convex setting.

Background and related work

Using a more tractable surrogate in the place of the actual target loss is a very old idea in statistics and machine learning, and some type of surrogate is used in almost any application. For instance, convex loss functions like the hinge or logistic loss are often used as surrogates for the discontinuous 0-1 classification loss. This avoids the computational intractability of minimizing the 0-1 loss (see, e.g., [START_REF] Arora | The hardness of approximate optima in lattices, codes, and systems of linear equations[END_REF], but note that this does not per se mitigate other difficulties with accessing the objective in the sense that we have discussed, e.g., limited access to training samples. Another classic example is empirical risk minimization [START_REF] Vapnik | The uniform convergence of frequencies of the appearance of events to their probabilities[END_REF] where the training loss is used as a surrogate for the real objective, the population loss. This more closely aligns with our motivation since it replaces the hard-to-evaluate population loss-which is defined by an unknown data distribution-with the training loss, which we can directly evaluate, compute gradients of, etc.

However, both of these examples rely on i.i.d. samples from the target distribution while our approach can exploit additional "side-information"-synthetic data, simulators, etc.to complement a smaller collection of i.i.d. samples. Moreover, the way in which surrogates are classically used differs from our approach: typically, the surrogate is directly minimized and the result is taken as an estimate of the minimum of the target objective, whereas we use the proxy to facilitate optimization using stochastic gradients from the objective itself. As such, the classic approach relies upon minima of the surrogate approximating minima of the target, while our method only requires the proxy and objective to have similar second derivatives. This is not necessarily easier to achieve, but it is easier in certain cases (see, e.g., Section 3.1), and it is much more amenable to applications involving synthetic data, simulators, etc., where the proxy loss could have different minima due to idiosyncrasies of the synthetic data or simulator.

In more closely related work, [START_REF] Hendrikx | Statistically preconditioned accelerated gradient method for distributed optimization[END_REF] propose an algorithm, SPAG, for a strongly convex setting, which builds upon a previous method, DANE [START_REF] Shamir | Communicationefficient distributed optimization using an approximate newton-type method[END_REF], and which resembles an accelerated version of Algorithm 1. Their guarantees are analogous to ours up to acceleration, however both SPAG and DANE require exact proximal-point updates, while our algorithm and analysis allow these updates to be inexact. Motivated by the problem of minimizing over a subset of parameters, [START_REF] Parpas | A multilevel proximal gradient algorithm for a class of composite optimization problems[END_REF] assumed that the surrogate subproblem is defined over a different space, with projection into that space and back only introducing a multiplicative error. In other related work, [START_REF] Chayti | Optimization with access to auxiliary information[END_REF] propose several algorithms that use a surrogate for non-convex optimization, one of which resembles a momentum variant of Algorithm 1. Mairal (2013a;2015) studies the "majorizationminimization" meta-algorithm, where an upper bound on the objective is minimized in each iteration. Our method can be viewed as an instance of majorization-minimization with the surrogate being used to define the upper bound on the loss. Mairal analyzes this method under the same

Algorithm 1 PROXYPROX 1: Input: initialization w 0 ∈ R d , stepsize η > 0 2: for k = 0, 1, . . . , K -1 do 3: Sample g k such that E[g k | w k ] = ∇L(w k ) 4:
Set w k+1 ≈ arg min w ϕ k (w) where

ϕ k (w) := g k , w + D L(w; w k ) + 1 2η w -w k 2 .
δ-similarity condition that we use and prove similar convergence rates, but the key difference between this work and ours is that Mairal's methods require a guaranteed upper bound on the objective-or on components of the objective in the case of finite sum structured problems-and require that this upper bound be exactly minimized in each iteration, whereas we only require an in-expectation upper bound and allow for only approximate minimization. In this way, our setting is more appropriate for the stochastic optimization problems that arise in machine learning.

Finally, in the most closely related work, Mairal (2013b) study a stochastic majorization-minimization approach which is applicable to the stochastic optimization problems we are interested in. However, Mairal's approach requires a guaranteed upper bound of the loss function evaluated on each sample and requires that the δ-similarity bound hold almost surely for each sample. In contrast, we only require an in-expectation upper bound and that the δ-similarity hold in expectation, making it easier to find a suitable surrogate. Furthermore, our convergence guarantees have a better scaling with the number of iterations, K.

The algorithm and its analysis

We now present our algorithm and prove convergence guarantees in the convex and strongly convex settings. We will also argue that our method can provide substantially better guarantees than SGD when L is sufficiently similar to L.

Setting and notation

Throughout this paper, we will use u to denote the Euclidean norm of the vector u. We will assume the minimum value of L, which we denote L * , is realized, with w * denoting an arbitrary minimizer. We use B 2 := E w 0 -w * 2 as shorthand to denote the distance between the (possibly random) initialization, w 0 , and this minimizer.

A function f is µ-strongly convex if for all u, v f (v) ≥ f (u) + ∇f (u), v -u + µ 2 u -v 2 , (4)
where, in the event that f is not differentiable, ∇f (u) denotes here an arbitrary subgradient of f . When this holds for µ = 0, f is merely convex. A function f is H-smooth if it is differentiable and ∇f is H-Lipschitz continuous or, equivalently, for all u, v

|f (v) -f (u) -∇f (u), v -u | ≤ H 2 u -v 2 . (5)
We define the Bregman divergence with potential ψ as

D ψ (u; v) := ψ(u) -ψ(v) -∇ψ(v), u -v . (6)
While Bregman divergences are typically defined only for strictly convex potentials ψ, in this paper we allow ψ to be any differentiable function with D ψ simply defined as in (6).

The key property of Bregman divergences for our purposes is the three-point identity:1 for all u, v, w

D ψ (u; v) -D ψ (u; w) -D ψ (w; v) = ∇ψ(v) -∇ψ(w), w -u . (7)
Our algorithm relies upon two main assumptions; one which captures the similarity between the objective L and the proxy L, and one which captures the level of noise in the stochastic gradients g k . First, we assume Assumption 2.1. The function h := L -L is differentiable and ∇h is δ-Lipschitz continuous.

That is, we require that our proxy L has similar curvature to the objective L. In light of (5), Assumption 2.1 implies that

|D h (u; v)| = D L (u; v) -D L(u; v) ≤ δ 2 u -v 2 .
This justifies the definition of ϕ k in Algorithm 1 because, when δ is small

ϕ k (w) ≈ g k , w + D L (w; w k ) + 1 2η w -w k 2 ,
which corresponds to (inexact, stochastic) proximal-point iterations directly on the objective L, which is well-known to have good performance (see, e.g., [START_REF] Barré | Principled analyses and design of first-order methods with inexact proximal operators[END_REF].

We note that Assumption 2.1 can accommodate nondifferentiable objectives, e.g., if L(w) = (w) + ψ(w) and L(w) = ˆ (w) + ψ(w) with and ˆ differentiable, because ψ does not appear in L-L. In this way, our algorithm is compatible, e.g., with non-smooth regularizers such as an L 1 penalty as long as this is incorporated into the proxy.

In addition to Assumption 2.1, we require the following:

Assumption 2.2. For each k, E[g k | w k ] = ∇L(w k ) and E g k -∇L(w k ) 2 w k ≤ σ 2 . If L is not differentiable but convex
, this holds for some subgradient in ∂L(w k ).

The assumption that the stochastic gradients g k are unbiased and have bounded variance is very common in the optimization literature [START_REF] Nemirovskii | Problem Complexity and Method Efficiency in Optimization[END_REF][START_REF] Dekel | Optimal distributed online prediction using mini-batches[END_REF][START_REF] Bubeck | Convex optimization: Algorithms and complexity[END_REF].

Convergence guarantees and proof sketch

We begin with our main result:

Theorem 2.3. Under Assumptions 2.1 and 2.2, let L be µ-strongly convex and let η ≤ 1 4δ . If for each k

E ∇ϕ k (w k+1 ) 2 ≤ µ 4η E w k+1 -w k 2 + G 2 ,
for some G, then a geometrically weighted average of the iterates of Algorithm 1, wK (see (25)), satisfies

EL( wK ) -L * ≤ 5B 2 8η 1 + 2ηµ 5 1-K + 2ησ 2 + G 2 µ ,
where the expectation is taken over the randomness in both the stochastic gradients and the selection of the iterates.

A detailed proof can be found in Appendices A and B. We start with a straightforward algebraic manipulation of ∇ϕ k .

Recalling that h = L -L, we write2

∇ϕ k (w k+1 ) = -∇h(w k+1 ) + ∇h(w k ) + ∇L(w k+1 ) + g k -∇L(w k ) + w k+1 -w k η . (8)
Next, (6) implies

L(w k+1 ) -L * = ∇L(w k+1 ), w k+1 -w * -D L (w * ; w k+1 ). (9)
Rearranging (8) allows us to express ∇L(w k+1 ) in terms of ∇ϕ k (w k+1 ) and the other quantities, which we substitute in the place of ∇L(w k+1 ) in ( 9). After several straightforward algebraic manipulations, including using the three-point identity (7), we arrive at the key identity

L(w k+1 ) -L * (10) = 1 2η w k -w * 2 -D h (w * ; w k ) - 1 2η w k+1 -w * 2 + D h (w * ; w k+1 )        (Telescope) + ∇L(w k ) -g k , w k+1 -w * (Noise) + ∇ϕ k (w k+1 ), w k+1 -w * (Inexact) + D h (w k+1 ; w k ) - 1 2η w k+1 -w k 2 -D L (w * ; w k+1 ) . (Cancel)
Although the expression looks unwieldy at first, it leads us directly towards a strategy for analyzing Algorithm 1.

The Telescope terms resemble those typically found in the analysis of convex optimization algorithms. In the course of proving Theorem 2.3, we, roughly speaking, average (10) over k, causing all but two of these terms to cancel out, and we show that under Assumption 2.1 with η ≤ 1 4δ , their total contribution to the error is small. The Noise term captures the effect of using noisy gradients g k rather than ∇L(w k ). Under Assumption 2.2, its expectation is at most 2ησ 2 + 1 4η E w k+1 -w k 2 . Importantly, since g k is unbiased, we can bound this in terms of w k+1 -w k rather than w k+1 -w * .

The Inexact term captures the effect of only approximately computing the proximal-point update. If we could just set w k+1 = arg min w ϕ k (w), the Inexact term would be zero, but this is impossible in practice. Instead, a simple application of Young's inequality bounds this term by

1 µ ∇ϕ k (w k+1 ) 2 + µ 4 w k+1 -w * 2 .
Nowhere above did we use any property of L besides those implied by Assumptions 2.1 and 2.2, but to address the Cancel terms, which we use to eliminate unwanted quantities introduced above, we finally use that in the µ-strongly convex setting, D L (w * ; w k+1 ) ≥ µ 2 w k+1 -w * 2 . So, when w k+1 is chosen so that ∇ϕ k (w k+1 ) 2 is sufficiently small, this allows us to cancel out everything in Noise and Inexact except for 2ησ 2 , and this leaves an extra µ 4 w k+1 -w * 2 to be incorporated into the Telescope terms to generate geometric shrinkage. Here, the strong convexity of L is crucially needed to cancel the dependence of the Inexact term on w k+1 -w * -without strong convexity, it is not immediately clear how the Inexact term can be controlled without setting ∇ϕ k (w k+1 ) = 0. Filling in the details and using standard arguments completes the proof. While Theorem 2.3 requires L to be µ-strongly convex, a similar guarantee can be achieved when L is convex via regularization: Theorem 2.4. Under Assumptions 2.1 and 2.2, let L be convex and let η ≤ 1 4δ . Then applying Algorithm 1 to the regularized objectives

L (µ) (w) := L(w) + µ 2 w -w 0 2 and L(µ) (w) := L(w) + µ 2 w -w 0 2 with µ = 1 ηK and any η ≤ 1 4δ ensures that if for each k E ∇ϕ k (w k+1 ) 2 ≤ 1 4η 2 K E w k+1 -w k 2 + G 2 for some G, then EL 1 K K k=1 w k -L * ≤ 9B 2 8ηK + 2ησ 2 + ηKG 2 .
The proof, which we defer to Appendix B, simply observes that L (µ) and L(µ) satisfy the conditions needed for Theorem 2.3, with larger µ leading to faster convergence. At the same time, when µ is small enough, optimizing L (µ) is essentially equivalent to optimizing L itself. We conclude by choosing µ to balance between these competing goals.

Finally, choosing η optimally in the strongly convex and convex settings yields the following complexity guarantee, which we prove in Appendix C:

Corollary 2.5. In the setting of Theorem 2.3, there is a universal constant c s.t. for any ε > 0, when G 2 ≤ 1 2 µε and

η = 5 µ(K -1) 1 + ln 5B 2 µ(K -1) ε ,
then the output, ŵ, of Algorithm 1 will have suboptimality at most EL( ŵ) -L * ≤ ε whenever

K ≥ c • 1 + δ µ + σ 2 µε ln e + (µ + δ)B 2 ε + σ 2 B 2 ε 2 .
Likewise, in the setting of Theorem 2.4, there is a universal constant c s.t. for any ε > 0, when

G 2 ≤ σ 2 K + δ 2 B 2 K 2 and η = min 1 4δ , B σ √ K ,
then the output, ŵ, of Algorithm 1 will have suboptimality at most EL( ŵ) -L * ≤ ε whenever

K ≥ c • δB 2 ε + σ 2 B 2 ε 2 .

Making ∇ϕ k small

Theorems 2.3 and 2.4 require that each update satisfy

E ∇ϕ k (w k+1 ) 2 ≤ µ 4η E w k+1 -w k 2 + G 2 ,
with µ = 1 ηK in the convex setting. This is a convenient inexactness criterion, because it consists entirely of quantities that can be evaluated at the time of execution, and similar bounds on the inaccuracy have been studied in the context of the classic proximal-point algorithm (see, e.g., [START_REF] Rockafellar | Augmented Lagrangians and applications of the proximal point algorithm in convex programming[END_REF][START_REF] Barré | Principled analyses and design of first-order methods with inexact proximal operators[END_REF]. Nevertheless, this raises the question of how difficult it is to find such a w k+1 . It turns out that this can be quite easy due to the 1 η -strong convexity of the subproblem objective ϕ k . For example, we could use stochastic gradient descent: Proposition 2.6. Let L be convex and H-smooth and let stochastic gradients be available for L with variance at most ρ 2 . Then for constants c, c , the output, ŵ, of

T ≥ c(1 + Hη) max ln c (1 + Hη) 2 µη , ρ 2 G 2 steps of SGD on ϕ k initialized at w k satisfies E ∇ϕ k ( ŵ) 2 ≤ µ 4η ŵ -w k 2 + G 2 .
The proof, which we defer to Appendix E, is based on the observation that L being H-smooth implies that ϕ k is

(H+ 1 η )-smooth, so ∇ϕ k ( ŵ) 2 ≤ 2(H+ 1 η )(ϕ k ( ŵ)-ϕ * k )
, and plugging in the standard suboptimality guarantee for T steps of SGD on a smooth and strongly convex objective eventually leads to the desired bound.

There is nothing particularly special about using SGD to compute the update w k+1 . An advantage of our method is that any algorithm for approximately minimizing ϕ k will suffice when L is smooth, so we can exploit structure in ϕ k by using a specialized optimization algorithm. There are also other, more exotic SGD variants that directly target a solution with a small gradient norm [START_REF] Allen-Zhu | How to make the gradients small stochastically: Even faster convex and nonconvex SGD[END_REF][START_REF] Foster | The complexity of making the gradient small in stochastic convex optimization[END_REF], and can have a better dependence on ρ 2 .

Comparison against directly optimizing L

To get a better sense of what Theorems 2.3 and 2.4 mean, we will compare them against what we could do by just directly optimizing the objective L using stochastic gradients. Suppose that Assumption 2.2 holds, the objective L is H-smooth, w 0 -w * 2 = B 2 , and L is either convex or µ-strongly convex. In this case, a natural and popular approach to minimizing L using a total of K stochastic gradients would be to execute K steps of SGD. Ignoring constant factors, it is well-known that this would yield an ε-accurate solution when [START_REF] Nemirovskii | Problem Complexity and Method Efficiency in Optimization[END_REF])

K ≥ HB 2 ε + σ 2 B 2 ε 2 (11) K ≥ H µ ln HB 2 ε + σ 2 µε . (12)
Instead, given a proxy L satisfying Assumption 2.1, we could use the same number of stochastic gradients from L to execute K steps of Algorithm 1. Ignoring constants and a logarithmic factor in the strongly convex setting, Corollary 2.5 shows that this yields an ε-accurate solution when

K ≥ δB 2 ε + σ 2 B 2 ε 2 (13) K ≥ δ µ ln δB 2 ε + σ 2 µε . (14)
Comparing ( 11) and ( 12) against ( 13) and ( 14), we see that Algorithm 1 essentially replaces SGD's dependence on H, the smoothness constant of L, with δ, the smoothness constant of h = L -L. So if a proxy L can be found with δ H, Algorithm 1 may have a much better guarantee than SGD.

On the other hand, the second, "statistical terms" of our guarantees match those of SGD. These statistical terms capture information-theoretic limits to optimizing L using only K stochastic gradients with variance σ 2 , and these are known to be unimprovable in the worst case [START_REF] Nemirovskii | Problem Complexity and Method Efficiency in Optimization[END_REF], so our algorithm is statistically optimal (up to a logarithmic factor in the strongly convex setting). In fact, long-standing lower bounds show that obtaining an ε-accurate solution using any algorithm that only uses K stochastic gradients from L requires at least [START_REF] Nemirovskii | Problem Complexity and Method Efficiency in Optimization[END_REF])

K ≥ HB 2 √ ε + σ 2 B 2 ε 2 K ≥ H µ ln HB 2 ε + σ 2 µε .
Therefore, when δ ≤ H √ ε in the convex setting or δ √ Hµ in the strongly convex case, then Algorithm 1 can surpass this lower bound by exploiting the proxy L in order to achieve a better guarantee than could be achieved by any algorithm that only uses K stochastic gradients from L.

Although Algorithm 1 only requires one stochastic gradient from L for each iteration, each update does require finding an approximate stationary point of the function ϕ k , which is more computationally expensive than one SGD update. That said, in the µ-strongly convex case (the conclusion is similar in the merely convex setting), if we ignore all constant and logarithmic factors, Corollary 2.5 shows that for K := δ µ + σ 2 µ , G 2 = µ , and η = 1 µK , Algorithm 1 will find an O( )-approximate minimizer of L. Furthermore, Proposition 2.6 (with ρ = σ) shows that the proximal subproblem in each iteration of Algorithm 1 can be solved to sufficient accuracy using T := K µK + Hσ 2 µ 2 K steps of SGD on ϕ k . Therefore, the total number of gradient update computations needed to implement Algorithm 1 is

T K = H µ + H µ • σ 2 µ .
In contrast, (12) implies that SGD directly on L would require S := H µ + σ 2 µ gradient update computations to reach an O( )-approximate minimizer of L.

So, SGD directly on L could be computationally cheaper than Algorithm 1 in terms of the total number of floating point operations needed since S < T K. However, crucially, Algorithm 1 only requires K stochastic gradients from L itself (the other gradients come from L), while SGD needs S > K. So to summarize, executing Algorithm 1 could require more computation than SGD on L, but it can make up for this by using fewer stochastic gradients from L, which could yield significant savings when L is harder to access.

Connection to Preconditioning

A useful way to understand Algorithm 1 is to consider the special case where the proxy L(w) = 1 2 w Pw + b w is quadratic. In this case, exactly solving the proximal subproblems in Algorithm 1 amounts to solving

g k + ∇ L(w k+1 ) -∇ L(w k ) + 1 η (w k+1 -w k ) = 0 ⇐⇒ w k+1 = w k -η(I + ηP) -1 g k . ( 15
)
That is, each update in Algorithm 1 amounts to a single preconditioned SGD update with preconditioning matrix (I + ηP) -1 . When a quadratic proxy that satisfies Assumption 2.1 can be found, then these updates with a fixed preconditioning matrix will be very effective at minimizing the objective.

However, for many objectives, the Hessian ∇ 2 L is not constant (nor close to constant) and so it cannot be wellapproximated everywhere by a matrix P and thus any fixed preconditioner may not be effective. In this more general case, we can view each step of Algorithm 1 as approximating a locally preconditioned SGD update with the preconditioning matrix

(I + η∇ 2 L(w k )) -1 . This is be- cause for sufficiently smooth L, ∇ L(w k+1 ) -∇ L(w k ) ≈ ∇ 2 L(w k )(w k+1 -w k ) so g k + ∇ L(w k+1 ) -∇ L(w k ) + 1 η (w k+1 -w k ) = 0 =⇒ w k+1 ≈ w k -η(I + η∇ 2 L(w k )) -1 g k . ( 16
)
Of course, this is not a perfect correspondence, but it demonstrates the role of L, which is effectively to adapt each update along the g k direction to the local curvature of L. Furthermore, when the curvature of L is close enough to that of L, this also indicates why Algorithm 1's updates should be more effective than SGD directly on L.

Applications

Algorithm 1 allows us to reduce the number of stochastic gradients needed from the objective L by approximately solving a sequence of strongly convex subproblems based on L. Accordingly, this approach is most advantageous when L is a faithful approximation of L in the sense that δ is small, and also when the model L is more accessible than the true objective L. This can arise in numerous situations, and in this section we will discuss several examples.

Regression and classification with costly labels

In supervised learning, we want to train using input-label pairs (x, y) ∼ D. Often, we can easily collect a large sample of i.i.d. inputs x 1 , x 2 , . . . , but obtaining labels for these inputs can be much more costly. For example, labelling whether an MRI scan contains a tumor may require consulting with (and likely paying) a trained medical technician. So, we would like to learn using as few labels as possible, and our algorithm can provide a means of doing this.

For inputs x ∈ R d and labels y ∈ R, least squares regression, where L(w) = 1 2 E (x,y)∼D [(w x -y) 2 ], is one of the main workhorses of statistical analysis. Conveniently, the Hessian of this loss, ∇ 2 L(w) = E x∼Dx [xx ], does not depend on the labels y, so we can construct a proxy L(w) = 1 2 w E x∼Dx [xx ]w that satisfies Assumption 2.1 with δ = 0 without using any labels! Thus, we can implement Algorithm 1 using just one labelled sample (x k , y k ) per iteration to estimate g k = (w x k -y k )x k , and then approximately minimizing ϕ k only requires x samples.

Similarly, for inputs x ∈ R d and binary labels y ∈ {0, 1}, the logistic regression loss function is

L(w) = -E (x,y)∼D y ln s(w x) + (1 -y) ln(1 -s(w x)) ,
where s(z) = 1/(1+exp(-z)) denotes the logistic function. As in the previous example, the Hessian of this loss is

∇ 2 L(w) = E (x,y)∼D s(w x)(1 -s(w x))xx ,
which again does not depend on the labels y. Therefore, we can construct a proxy loss by simply assigning all points the label y = 1 (or assigning them any arbitrary labels):

L(w) := -E x∼Dx [ln s(w x n )].
This ensures ∇ 2 L = ∇ 2 L so that Assumption 2.1 holds with δ = 0, despite using no labels. So again, Algorithm 1 can be implemented using as little as one label per iteration, plus sufficient x samples to approximately minimize ϕ k .

Synthetic data, simulators, and meta-learning

In many machine learning problems of interest, it is difficult to collect i.i.d. samples from the target distribution. For instance, as discussed in the previous section, obtaining a label in medical imaging applications may require a consultation with a radiologist and perhaps even an additional biopsy in difficult cases. This can make it costly and timeconsuming to obtain a large training set. To mitigate this issue, there has been great interest in using synthetic data as a complement to or substitute for costly samples from the target distribution. Recent advances in generative models like generative adversarial networks (GANs), variational autoencoders, diffusion models, etc. have raised the possibility of using a small seed of real data to generate an essentially unlimited quantity of synthetic training data. For instance, Frid-Adar et al. ( 2018) used a GAN to augment a small dataset of 182 CT scans of liver lesions, which allowed them to train a classifier with better performance. Similarly, meta-learning [START_REF] Hospedales | Meta-learning in neural networks: A survey[END_REF], multitask learning [START_REF] Zhang | A survey on multi-task learning[END_REF], and transfer learning [START_REF] Weiss | A survey of transfer learning[END_REF] are all efforts to improve performance on a given task by leveraging data collected for other, related tasks. For example, when learning to classify images of giraffes, it is common to use another image dataset like ImageNet for pretraining or joint training-even though ImageNet does not contain labelled images of giraffes.

In addition to possible difficulties in obtaining the training data, evaluating the loss or its gradient can be expensive. In robotics, evaluating the reward achieved by a given policy might require waiting for a physical robot to execute the policy for a period of time, and waiting seconds or minutes for each function or gradient evaluation can greatly slow training. In addition, training something like a self-driving car in the real world carries the risk of damaging the vehicle or injuring pedestrians. As a result, it is very common to use simulators like MuJoCo [START_REF] Todorov | A physics engine for model-based control[END_REF] and the OpenAI Gym (Brockman et al., 2016), which make it possible to train models in a virtual environment, allowing for faster, risk-free processing.

Synthetic data and simulators are often used to define a surrogate loss function that is directly minimized to train the model. In meta-learning, multitask learning, or transfer learning, other tasks are usually used as a prior that biases training towards models that also perform well on all tasks simultaneously. Therefore, all of these approaches will only work to the extent that the minima of the surrogate losses correspond to minima of the target loss L itself. Indeed, recent work has cast some doubt on the utility of, e.g., ImageNet pretraining [START_REF] He | Rethinking imagenet pretraining[END_REF][START_REF] Kornblith | Do better imagenet models transfer better[END_REF] and robots trained on simulators frequently generalize poorly to the real world [START_REF] Zhao | Sim-to-real transfer in deep reinforcement learning for robotics: a survey[END_REF], presumably because these surrogates can bias the model towards worse solutions. In contrast, our approach does not rely on minima of the proxy having any particular correspondence with those of the target; instead, we use the curvature of the proxy to speed up optimization of the objective itself.

Learning with public and private data

Another natural setting where our algorithm could be useful is when the objective L is based on private data. One approach to learning from this data while maintaining rigorous privacy guarantees is to use a differentially private algorithm to generate a synthetic dataset, which could then be used freely as described in the previous section (see, e.g., [START_REF] Amin | Differentially private covariance estimation[END_REF][START_REF] Bowen | Comparative study of differentially private synthetic data algorithms from the NIST PSCR differential privacy synthetic data challenge[END_REF] Our method might also be useful when a stock of non-private data is be available in addition to the sensitive data. For example, some patients could elect to provide their medical data to researchers without restriction, while others only do so on the condition that their privacy is rigorously protected. Similarly, different individuals may be subject to different regulatory regimes-data from E.U. residents must be used in accordance with more stringent GDPR rules that do not apply to, e.g., US residents. In these cases, one could use just the non-private data to learn, and avoid the trouble of using the protected data, but this could be suboptimal if a large fraction of the total data is private. Instead, we could use the non-private data to define L while using a privacypreserving mechanism to compute stochastic gradients of L from the protected samples [START_REF] Dwork | The algorithmic foundations of differential privacy[END_REF]. When the private and non-private samples have similar distributions, we would expect Assumption 2.1 to hold with small δ, meaning that relatively few iterations of Algorithm 1 will be needed to achieve low error, and therefore the private data would need to be accessed fewer times. Since privacy guarantees degrade with the number of accesses, this could allow for better performance and better privacy guarantees than would be possible using the private data alone.

Handling distribution shifts

Often, we would like to learn in the face of distribution shift-a situation where the training distribution differs from the test distribution. For example, the distribution of inputs and outputs could be changing over time, and the training set is collected today and then the model is deployed tomorrow. In recent work, [START_REF] Lazaridou | Mind the gap: Assessing temporal generalization in neural language models[END_REF] argue that the performance of neural language models can degrade over time for this very reason. To keep a model up-to-date, one could periodically collect a new dataset and train a whole new model. However, collecting a large enough quantity of high-quality data could be expensive and time consuming, and it would be better to take advantage of the older data that is already available. To apply our algorithm, we could have L be the loss on the current data distribution, while out-of-date data could be used to define L. As long as the distribution shift is not too dramatic, we could expect Assumption 2.1 to hold with a relatively small δ, and so to implement our method, only a small amount of new data would needed.

Distribution shifts can also arise when the training distribution is stationary but biased. For example, survey data is often affected by selection, acquiescence, or social desirability biases which can lead to differences between survey results and reality. Relatedly, when learning from a data source in which one group or class is rare, we may upweight the loss on the rare group in order to prevent the model from learning to ignore it, which corresponds to targeting a test distribution where that group's prevalence is higher. In either case, our algorithm could be applied as described above using relatively few unbiased samples.

Experiments

To exhibit the effectiveness of our method in practice, we conducted several experiments on realistic problems.

Logistic regression: First, we consider a simple binary logistic regression problem with the 'mushrooms' dataset, which consists of 8124 samples of dimension 112. In this experiment, the objective L is the logistic regression loss evaluated on the training data plus an 2 regularizer with weight µ = 10 -6 H, where H is the smoothness parameter of the objective. Stochastic gradients for L are calculated by evaluating the gradient on a minibatches of size 256 or 1024 drawn uniformly with replacement. To simulate a situation where labels are hard to come by as discussed in Section 3.1, we define the proxy L by assigning random labels to the input features. The results, depicted in Figure 1, show that our method PROXYPROX has competitive performance compared with performing SGD directly on the objective. We are also competitive with the AUXMOM and AUXMVR algorithms proposed by [START_REF] Chayti | Optimization with access to auxiliary information[END_REF], which are momentum and variance-reduced analogues of PROXYPROX. In Figure 1, the x-axis counts how many stochastic gradients from L have been used by each method, so this comparison holds constant the amount of access to L, although the computational cost of PROXYPROX, AUXMOM and AUXMVR is higher due to the approximate proximal-point updates in each iteration. We approximately solve each of the ϕ k subproblems using gradient descent.

For each algorithm in the experiment, we selected hyperparameters like η using grid search to minimize the loss after 250 steps and ran the best option for the full 1000 iterations. We set AUXMOM and AUXMVR's a = 0.1, which corresponds to the standard momentum setting of 0.9.

ResNet-18: In this experiment, summarized in Figure 2, we train a ResNet-18 network on CIFAR-10, defining L using a 2560 subset of the images, and using the full train dataset of 50000 images for L. While L corresponds to the training loss, we also plot the test loss to show that our method does not overfit. We compare running our method with either 20 and 40 iterations of SGD (equivalent to 1 epoch and 2 epochs on the L data) to approximately minimize ϕ k , with the stepsize tuned by grid search and ultimately set to 0.01.

Figure 2 indicates that a single pass over the subsampled data seems sufficient to minimize ϕ k . We use minibatch size of 128 for both g k and the SGD updates to minimize ϕ k . Our method can be better than ADAMW for the first ∼ 25 epochs, and it improves more convincingly over SGD. We use standard stepsizes: 0.1 for SGD, 0.001 for ADAMW, and weight decay of 0.1 for ADAMW, which gave the best test accuracy in a grid search, and all methods also used cosine annealing. We note that just training on L using SGD plateaus at test accuracy 65 after just a few passes, so, our method's improvement is indeed coming from using L jointly with L.

Open problems and future directions

In the convex setting, we show in Section 2.2 that under Assumptions 2.1 and 2.2, Algorithm 1 can take advantage of the proxy to achieve a convergence rate comparable to SGD on a δ-smooth function-despite the fact that L may not be smooth, and even if it is smooth, its smoothness parameter may be much larger than δ. However, the complexity guarantees in Corollary 2.5 require knowledge of problem parameters, including δ which may be hard to estimate in practice. In future work, it may therefore be useful to explore whether similar performance could be achieved using a method that is adaptive to unknown problem parameters. Also, the uniform bound on the stochastic gradient variance in Assumption 2.2 can be unrealistic, and in other contexts it can be possible to show similar results under the weaker assumption that the gradient variance is bounded only at w * [START_REF] Bach | Non-asymptotic analysis of stochastic approximation algorithms for machine learning[END_REF][START_REF] Dragomir | Fast stochastic Bregman gradient methods: Sharp analysis and variance reduction[END_REF].

Going beyond the convex setting is also important in machine learning, where many models used in practice, like neural networks, give rise to non-convex training objectives. Conceptually, the idea of using the proxy L to guide optimization of the objective L seems sound even when these functions are non-convex, and our CIFAR-10 experiment shows that it can work well. However, there is still work to do in proving theoretical guarantees. We can show that Algorithm 1 converges to the vicinity of a stationary point: Theorem 5.1. Under Assumptions 2.1 and 2.2, with L differentiable but potentially non-convex, let η ≤ 1 4δ . If

E ϕ k (w k+1 ) 2 ≤ 7 16η 2 E w k+1 -w k 2 + 1 8 ∇L(w k+1 ) 2 ,
and Eϕ k (w k+1 ) ≤ Eϕ k (w k ) for each k, then

1 K K k=1 E ∇L(w k ) 2 ≤ 48(L(w 0 ) -L * ) ηK + 8σ 2 .
We prove this in Appendix D. Taking η = 1 4δ , this first term becomes O(δ(L(w 0 ) -L * )/K), which matches the rate of gradient descent on a δ-smooth objective, paralleling our results in the convex setting. However, the noise term remains Ω(σ 2 ) regardless of η or K. It is not clear whether this is merely a shortcoming in our analysis, or if this actually reflects the worst case performance of Algorithm 1. To our knowledge, a satisfying analysis of stochastic proximalpoint methods in the non-convex setting is lacking more generally, so this situation may not be unique to our method.

A. Proof of Theorem 2.3

Our analysis requires that h = L -L be differentiable under Assumption 2.1, but when L is convex, it is not necessary for L itself to be differentiable. In particular, Assumption 2.2 requires that the expectation of each stochastic gradient g k must be a subgradient of L at w k . Accordingly, in the following proofs, we will use ∇L(w k ) := E[g k | w k ] as notation to denote the particular subgradient corresponding to the expectation of g k . Since h = L -L is differentiable, we can use this, in turn, to define ∇ L(w k ) := ∇L(w k ) -∇h(w k ). Finally, since 0 ∈ ∂L(w * ), we denote ∇L(w * ) := 0 even when L is not differentiable at its minimizer.

Lemma A.1. Let w * ∈ arg min w L(w) and η ≤ 1 4δ . Then under Assumptions 2.1 and 2.2, for any w k+1 , w k

E[L(w k+1 ) -L * ] ≤ 1 2η E w k -w * 2 -ED h (w * ; w k ) - 1 2η E w k+1 -w * 2 + ED h (w * ; w k+1 ) + E ∇ϕ k (w k+1 ), w k+1 -w * + 2ησ 2 - 1 4η E w k+1 -w k 2 -ED L (w * ; w k+1 ),
where the expectation is taken over the randomness in g k .

Proof. Computing ∇ϕ k and rearranging yields

∇ϕ k (w k+1 ) = ∇ L(w k+1 ) + g k -∇ L(w k ) + 1 η (w k+1 -w k ) =⇒ ∇L(w k+1 ) = ∇h(w k+1 ) -∇h(w k ) + ∇ϕ k (w k+1 ) + ∇L(w k ) -g k - 1 η (w k+1 -w k ). (17)
Also, (6) immediately implies L(w k+1 ) -L * = ∇L(w k+1 ), w k+1 -w * -D L (w * ; w k+1 ).

After substituting ( 17) into (18), we derive

L(w k+1 ) -L * = ∇L(w k+1 ), w k+1 -w * -D L (w * ; w k+1 ) = ∇h(w k+1 ) -∇h(w k ) + ∇ϕ k (w k+1 ) + ∇L(w k ) -g k + 1 η (w k -w k+1 ), w k+1 -w * -D L (w * ; w k+1 ) = D h (w * ; w k+1 ) -D h (w * ; w k ) + D h (w k+1 ; w k ) + ∇ϕ k (w k+1 ) + ∇L(w k ) -g k , w k+1 -w * + 1 2η w k -w * 2 - 1 2η w k+1 -w * 2 - 1 2η w k+1 -w k 2 -D L (w * ; w k+1 ), (19)
where the third equality uses the three-point property with D h , (7). To proceed, Assumption 2.2, Young's inequality, and η ≤ 1 4δ , together imply

E ∇L(w k ) -g k , w k+1 -w * = E ∇L(w k ) -g k , w k+1 -w k ≤ ηE ∇L(w k ) -g k 2 1 -2ηδ + 1 -2ηδ 4η w k+1 -w k 2 ≤ 2ησ 2 + 1 -2ηδ 4η w k+1 -w k 2 ,
where the expectation is taken over the randomness in g k . Since h is δ-smooth, in light of (5), |D h (w;

w k )| ≤ δ 2 w -w k 2 , so E[D h (w k+1 ; w k ) + ∇L(w k ) -g k , w k+1 -w * ] ≤ 2ησ 2 + 1 4η E w k+1 -w k 2 .
Taking the expectation of ( 19) and plugging this in completes the proof.

Theorem 2.3. Under Assumptions 2.1 and 2.2, let L be µ-strongly convex and let η ≤ 1 4δ . If for each k

E ∇ϕ k (w k+1 ) 2 ≤ µ 4η E w k+1 -w k 2 + G 2 ,
for some G, then a geometrically weighted average of the iterates of Algorithm 1, wK (see (25)), satisfies

EL( wK ) -L * ≤ 5B 2 8η 1 + 2ηµ 5 1-K + 2ησ 2 + G 2 µ ,
where the expectation is taken over the randomness in both the stochastic gradients and the selection of the iterates.

Proof. For brevity, denote

∆ k := 1 2η w k -w * 2 -D h (w * ; w k ).
Because L is strongly convex, it follows immediately from ( 4) and ( 6) that D L (w * ; w k+1 ) ≥ µ 2 w k+1 -w * 2 . Also, in light of (5), it follows under Assumption 2.1 that

|D h (w * ; w k+1 )| ≤ δ 2 w k+1 -w * 2 ≤ 1 8η w k+1 -w * 2 .
This implies that

0 ≤ ∆ k+1 ≤ 5 8η k+1 -w * 2 ≤ 5 4ηµ D L (w * ; w k+1 ).
Therefore, by Lemma A.1 and Young's inequality

E[L(w k+1 ) -L * ] ≤ E∆ k -E∆ k+1 + 2ησ 2 + E ∇ϕ k (w k+1 ), w k+1 -w * - 1 4η E w k+1 -w k 2 -ED L (w * ; w k+1 ) ≤ E∆ k -1 + 2ηµ 5 E∆ k+1 + 2ησ 2 + E ∇ϕ k (w k+1 ), w k+1 -w * - 1 4η E w k+1 -w k 2 - µ 4 E w k+1 -w * 2 . ( 20
)
So by Young's inequality and the assumed upper bound

∇ϕ k (w k+1 ) 2 ≤ µ 4η E w k+1 -w k 2 + G 2 (21)
we conclude

E ∇ϕ k (w k+1 ), w k+1 -w * ≤ 1 µ E ∇ϕ k (w k+1 ) + µ 4 E w k+1 -w * 2 (22) ≤ 1 4η E w k+1 -w k 2 + µ 4 E w k+1 -w * 2 + G 2 µ . (23)
Plugging this into (20) yields

E[L(w k+1 ) -L * ] ≤ E∆ k -1 + 2ηµ 5 E∆ k+1 + 2ησ 2 + G 2 µ . (24)
As is typical in strongly convex stochastic optimization, we now introduce weights α k = (1 + 2ηµ 5 ) k-1 (see Section 5 in [START_REF] Nesterov | Confidence level solutions for stochastic programming[END_REF], for a similar weighting) with sum A K := K k=1 α k and we will attempt to upper bound the suboptimality of the weighted iterate

wK := 1 A K K k=1 α k w k . (25)
Therefore, for z ≥ 0 ln

(1 + z) = z 0 g (t)dt ≥ z 0 (1 -t)dt = z - 1 2 z 2 .
Thus,

1 + 2ηµ 5 1-K = exp -(K -1) ln 1 + 2ηµ 5 ≤ exp - 2(K -1)ηµ 5 1 - ηµ 5
therefore, as long as we choose η ≤ 5 2µ , then it follows that

5B 2 8η 1 + 2ηµ 5 1-K ≤ 5B 2 µ(K -1) 4 exp - (K -1)ηµ 5 .
Therefore, if all of the following constraints on η are satisfied:

η ≤ 1 4δ (27) η ≤ 5 2µ (28) η ≥ 1 2µ(K -1) (29) η ≥ 5 µ(K -1) ln 5B 2 µ(K -1) 4ε (30)
then this first term in the guarantee from Theorem 2.3 is upper bounded by ε. We therefore set

η = 5 µ(K -1) 1 + ln 5B 2 µ(K -1) 4ε . (31)
This value always satisfies ( 29) and ( 30), and we will proceed to choose K large enough that it also satisfies ( 27) and ( 28).

In particular, we require K such that

K ≥ 1 + 2 max 1, 10δ µ 1 + ln 5B 2 µ(K -1) 4ε . ( 32
)
Plugging the stepsize (31) into the second term of the guarantee from Theorem 2.3 yields

2ησ 2 = 10σ 2 µ(K -1) 1 + ln 5B 2 µ(K -1) 4ε .
So, to make this smaller than ε, we require

K ≥ 1 + 10σ 2 µε 1 + ln 5B 2 µ(K -1) 4ε .
Combining this with (32), we conclude that with η set according to (31), if

K ≥ 10 max 1, δ µ , σ 2 µε 1 + ln 5B 2 µK 4ε (33)
and, in addition G ≤ √ µε, then the expected suboptimality is at most 3ε, so using ε = ε/3 instead above ensures accuracy ε.

To conclude the proof, we note for A, B > 0, the inequality

X ≥ A(1 + ln(BX))
is satisfied by setting X ≥ A(1 + 4 ln(e + AB)).

To show this, we first note that because z → √ z is a concave function on z > 0, we have the inequality

√ z ≤ √ 1 + (z -1) d dz √ z z=1 = 1 + z -1 2 = 1 + z 2 .
This implies that for z > 0 d dz

√ z = 1 2 √ z ≥ 1 1 + z = d dz ln(1 + z), (34)
Finally, since √ 0 = ln(1 + 0), this implies that ln(1 + z) ≤ √ z for all z ≥ 0.

So, suppose that X = αA(1 + 4 ln(e + AB)).

for some α ≥ 1. Then,

A(1 + ln(BX)) = A(1 + ln(α) + ln(AB) + ln(1 + 4 ln(e + AB)))

≤ A(1 + ln(α) + ln(e + AB) + 4 ln(e + AB))

≤ A(1 + ln(α) + 3 ln(e + AB))

≤ 1 α X + ln(c)A ≤ 1 + ln(α) α X ≤ X,
where we used that 1 + ln(α) ≤ α for all α ≥ 1.

Therefore, for a sufficiently large constant c, setting Rearranging and substituting ϕ k , this implies

K = c • 1 + δ µ +
L(w) -L(w k ) ≤ ϕ k (w) -ϕ k (w k ) - 1 -ηδ 2η w -w k 2 + ∇L(w k ) -g k , w -w k .
Next we use Young's inequality with Assumption 2.2 to upper bound

E ∇L(w k ) -g k , w -w k ≤ ηE ∇L(w k ) -g k 2 1 -2ηδ + 1 -2ηδ 4η E w -w k 2 ≤ ησ 2 1 -2ηδ + 1 -2ηδ 4η E w -w k 2 .
The fact that η ≤

≥ 1 4 E ∇L(w) 2 -E ∇ϕ k (w) 2 -δ 2 E w -w k 2 -σ 2 ≥ 1 8 E ∇L(w) 2 - 1 2η 2 E w -w k 2 -σ 2 ,
where we used the upper bound on E ∇ϕ k (w) 2 and the fact that η ≤ 1 4δ . Rearranging completes the proof. steps of SGD on ϕ k initialized at w k satisfies

E ∇ϕ k ( ŵ) 2 ≤ µ 4η ŵ -w k 2 + G 2 .
Proof. Since ϕ k (w) = L(w) + 1 2η ww k 2 up to an additional affine term, it is easy to see that ϕ k is 1 η -strongly convex and (H + 1 η )-smooth. Therefore, for any w

∇ϕ k (w) 2 ≤ 2 H + 1 η (ϕ k (w) -ϕ k (w * ϕ )), (40)
where w * ϕ = arg min w ϕ k (w). In addition, the output, ŵ, of SGD on ϕ k using the optimal stepsize and initialized at w k is, for universal constants c, c ≥ 1 [START_REF] Nemirovskii | Problem Complexity and Method Efficiency in Optimization[END_REF])

E[ϕ k ( ŵ) -ϕ k (w * ϕ )] ≤ c H + 1 η w k -w * ϕ 2 exp - T c (1 + Hη) + cηρ 2 T . (41)
So, by the relaxed triangle inequality and the strong convexity of ϕ k

w k -w * ϕ 2 ≤ 2 ŵ -w k 2 + 2 ŵ -w * ϕ 2 (42) ≤ 2 ŵ -w k 2 + 4η ϕ k ( ŵ) -ϕ k (w * ϕ ) (43)

Figure 1 .Figure 2 .

 12 Figure 1. Binary logistic regression on 'mushrooms'. Left: batch size 1024; middle: batch size 256; right: performance of PROXYPROX with L subsampled from L, stochastic gradients of L computed with batch size 256.

  1. Under Assumptions 2.1 and 2.2, let η ≤ 1 4δ . Then for any w such that Eϕ k(w) ≤ Eϕ k (w k ), E[L(w) -L(w k )] ≤ -1 4η E ww k 2 + 2ησ 2 . Proof. Under Assumption 2.1, h = L -L is δ-smooth. So, L(w) -L(w) ≤ L(w k ) -L(w k ) + δ 2 ww k 2 + ∇L(w k ) -∇ L(w k ), ww k .

Proof.

  By the relaxed triangle inequality, for any vectors a, b, c, d,a 2 = a + b + c + d -b -c -d 2 ≤ 4( a + b + c + d 2 + b 2 + c 2 + d 2 ).(36)Furthermore,∇ϕ k (w k+1 ) = ∇ L(w k+1 ) -∇L(w k+1 ) -∇ L(w k ) + ∇L(w k ) + g k -∇L(w k ) + ∇L(w k+1 ) + 1 η (w k+1 -w k ). (37)So (36) witha = ∇L(w k+1 ) b = -∇ϕ k (w k+1 ) c = ∇ L(w k+1 ) -∇L(w k+1 ) -∇ L(w k ) + ∇L(w k ) d = g k -∇L(w k ) so a + b + c + d = 1 η (w k -w k+1), implies that under Assumptions 2.1 and 2∇ϕ k (w) 2 -E g k -∇L(w k ) 2 -E ∇ L(w) -∇L(w) -∇ L(w k ) + ∇L(w k ) 2

  1 4δ and Eϕ k (w) ≤ Eϕ k (w k ) completes the proof. Lemma D.2. Under Assumptions 2.1 and 2.2, let L be differentiable and let η ≤ 1 4δ . If w k+1 satisfies E ϕ k (w k+1 ) 2 ≤

				7 16η 2 E w k+1 -w k	2 +	1 8	∇L(w k+1 ) 2
	then	-	1 4η	w k+1 -w k	2 ≤ -	η 48	E ∇L(w k+1 )	2 +	ησ 2 6	.

  Theorem 5.1. Under Assumptions 2.1 and 2.2, with L differentiable but potentially non-convex, let η ≤ 1 4δ . IfE ϕ k (w k+1 ) 2 ≤ and Eϕ k (w k+1 ) ≤ Eϕ k (w k ) for each k, then . Proof of Proposition 2.6Proposition 2.6. Let L be convex and H-smooth and let stochastic gradients be available for L with variance at most ρ 2 . Then for constants c, c , the output, ŵ, of T ≥ c(1 + Hη) max ln c (1 + Hη)

							µη	2	,	ρ 2 G 2
					7 16η 2 E w k+1 -w k	2 +	1 8	∇L(w k+1 ) 2 ,
	1 K	K k=1	E ∇L(w k ) 2 ≤	48(L(w 0 ) -L * ) ηK	+ 8σ 2 .
	Proof. By Lemmas D.1 and D.2
	E[L(w k+1 ) -L(w k )] ≤ -	1 4η	E w k+1 -w k	2 + 2ησ 2 ≤ -	η 48	E ∇L(w k+1 )	2 +	ησ 2 6	.	(38)
	Therefore, rearranging and averaging over k yields
	1 K										ηK	+ 8σ 2	(39)
	which completes the proof.

K k=1 E ∇L(w k ) 2 ≤ 48(L(w 0 ) -L(w K )) ηK + 8σ 2 ≤ 48(L(w 0 ) -L * ) E

Inria, Ecole Normale Supérieure, PSL Research University, Paris France

Samsung AI Center, Cambridge, UK. Work done while at CNRS, Ecole Normale Supérieure, Inria. Correspondence to: Blake Woodworth <blakewoodworth@gmail.com>, Konstantin Mishchenko <konsta.mish@gmail.com>, Francis Bach <francis.bach@inria.fr>.

This does not require convexity(Chen & Teboulle, 1993, Lemma 3.1).

Here, for simplicity, assume L and L are differentiable, but see Appendix A for the general case.

Acknowledgements

We acknowledge support from the French government under the management of the Agence Nationale de la Recherche as part of the "Investissements d'avenir" program, reference ANR-19-P3IA-0001 (PRAIRIE 3IA Institute), as well as from the European Research Council (grant SEQUOIA 724063).

By the convexity of L and (24),

Finally, we lower bound

, which completes the proof.

B. Proof of Theorem 2.4

Theorem 2.4. Under Assumptions 2.1 and 2.2, let L be convex and let η ≤ 1 4δ . Then applying Algorithm 1 to the regularized objectives L (µ) (w) := L(w) + µ 2 ww 0 2 and L(µ) (w) := L(w) + µ 2 ww 0 2 with µ = 1 ηK and any η ≤ 1 4δ ensures that if for each k

for some G, then

Proof. By construction, L (µ) is µ-strongly convex, and yet L (µ) -L(µ) = h, so this difference remains δ-smooth. Therefore, our analysis from Theorem 2.3 can be applied to L (µ) . Denote w * µ := arg min w L (µ) (w). In the course of proving Theorem 2.3, we showed in ( 24) that

where

Therefore, by the convexity of

Finally,

where for the last line we used the convexity of L. So,

When combined with (26), this means

and plugging in our choice of µ = 1 ηK completes the proof.

C. Proof of Corollary 2.5

Corollary 2.5. In the setting of Theorem 2.3, there is a universal constant c s.t. for any ε > 0, when G 2 ≤ 1 2 µε and

then the output, ŵ, of Algorithm 1 will have suboptimality at most EL( ŵ) -L * ≤ ε whenever

Likewise, in the setting of Theorem 2.4, there is a universal constant c s.t. for any ε > 0, when

then the output, ŵ, of Algorithm 1 will have suboptimality at most EL( ŵ) -L * ≤ ε whenever

Proof. The complexity guarantee in the convex case follows immediately after plugging the chosen η into the guarantee of Theorem 2.4 and then choosing K large enough that the expected suboptimality is less than ε.

Moving on to the strongly convex setting, to apply Theorem 2.3, we require η ≤ 1 4δ . Suppose for now that η also satisfies

.

Then the first term in the guarantee from Theorem 2.3 is at most

.

Next, we note that for g(z) = ln(1 + z), and any z ≥ 0

Finally, we note that our choice of T (51) satisfies the condition (45) because L being H smooth implies that L is (H + δ)-smooth and therefore