 # How Good is SGD with Random Shuffling?

We study the performance of stochastic gradient descent (SGD) on smooth and strongly-convex finite-sum optimization problems. In contrast to the majority of existing theoretical works, which assume that individual functions are sampled with replacement, we focus here on popular but poorly-understood heuristics, which involve going over random permutations of the individual functions. This setting has been investigated in several recent works, but the optimal error rates remains unclear. In this paper, we provide lower bounds on the expected optimization error with these heuristics (using SGD with any constant step size), which elucidate their advantages and disadvantages. In particular, we prove that after k passes over n individual functions, if the functions are re-shuffled after every pass, the best possible optimization error for SGD is at least Ω(1/(nk)^2+1/nk^3), which partially corresponds to recently derived upper bounds, and we conjecture to be tight. Moreover, if the functions are only shuffled once, then the lower bound increases to Ω(1/nk^2). Since there are strictly smaller upper bounds for random reshuffling, this proves an inherent performance gap between SGD with single shuffling and repeated shuffling. As a more minor contribution, we also provide a non-asymptotic Ω(1/k^2) lower bound (independent of n) for cyclic gradient descent, where no random shuffling takes place.

Comments

There are no comments yet.

## Authors

##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

We consider variants of stochastic gradient descent (SGD) for solving unconstrained finite-sum problems of the form

 minx∈XF(x) = 1nn∑i=1fi(x) , (1)

where is some Euclidean space (or more generally some real Hilbert space), is a strongly convex function, and each individual function

is smooth (with Lipschitz gradients) and Lipschitz on a bounded domain. Such problems are extremely common in machine learning applications, which often boil down to minimizing the average loss over

data points with respect to a class of predictors parameterized by a vector

. When is large, perhaps the most common approach to solve such problems is via stochastic gradient descent, which initializes at some point in and involves iterations of the form where . The majority of existing theoretical works assume that each is sampled independently across iterations (also known as with-replacement sampling). For example, if it is chosen independently and uniformly at random from , then , so the algorithm can be seen as a noisy version of exact gradient descent on (with iterations of the form ), which greatly facilitates its analysis.

However, this straightforward sampling approach suffers from practical drawbacks, such as requiring truly random data access and hence longer runtime. In practice, it is quite common to use without-replacement sampling heuristics, which utilize the individual functions in some random or even deterministic order (see for example [2, 3, 9, 13, 14, 1, 4]). Moreover, to get sufficiently high accuracy, it is common to perform several passes over the data, where each pass either uses the same order as the previous one, or some new random order. The different algorithm variants we study in this paper are presented as Algorithms 1 to 4 below. We assume that all algorithms take as input the functions , a step size parameter (which remains constant throughout the iterations), and an initialization point . The algorithms then perform

passes (which we will also refer to as epochs) over the individual functions, but differ in their sampling strategies:

• Algorithm 1 (SGD with random reshuffling) chooses a new permutation of the functions at the beginning of every epoch, and processes the individual functions in that order.

• Algorithm 2 (SGD with single shuffling) uses the same random permutation for all epochs.

• Algorithm 3 (Cyclic SGD) performs passes over the individual functions, each in the same fixed order (which we will assume without loss of generality to be the canonical order )

In contrast, Algorithm 4 presents SGD using with-replacement sampling, where each iteration an individual function is chosen uniformly and independently.

These without-replacement sampling heuristics are often easier and faster to implement in practice. In addition, when using random permutations, they often exhibit faster error decay than with-replacement SGD. A common intuitive explanation for this phenomenon is that random permutations force the algorithm to touch each individual function exactly once during each epoch, whereas with-replacement makes the algorithm touch each function once only in expectation. However, theoretically analyzing these sampling heuristics has proven to be very challenging, since the individual iterations are no longer statistically independent.

In the past few years, some progress has been made in this front, and we summarize the known results on the expected optimization error (or at least what these results imply111For example, some of these papers focus on bounding where is the minimum of , rather than the expected optimization error . However, for strongly convex and smooth functions, and are the same up to the strong convexity and smoothness parameters, see for example .), as well as our new results, in Table 1. First, we note that for SGD with replacement, classic results imply an optimization error of after stochastic iterations, and this is known to be tight (see for example ). For SGD with random reshuffling, better bounds have been shown in recent years, generally implying that when the number of epochs is sufficiently large, such sampling schemes are better than with-replacement sampling, with optimization error decaying quadratically rather than linearly in . However, the optimal dependencies on and other problem-dependent parameters remain unclear (HaoChen and Sra  point out that for , one cannot hope to achieve worst-case error smaller than , but for not much is known). Some other recent theoretical works on SGD with random reshuffling (but under somewhat different settings) include [13, 16]. For cyclic SGD, an upper bound was shown in , as well as a matching asymptotic lower bound in terms of . For SGD with single shuffling, we are actually not aware of a rigorous theoretical analysis. Thus, we only have the upper bound trivially implied by the analysis for cyclic SGD, and for , the upper bound implied by the analysis for random reshuffling (since in that case there is no distinction between single shuffling and random reshuffling). Indeed, for single shuffling, even different epochs are not statistically independent, which makes the analysis particularly challenging.

In this paper, we provide lower bounds on the expected optimization error of SGD with these sampling heuristics, which complement the existing upper bounds and provide further insights on the advantages and disadvantages of each. We focus on constant-step size SGD, as it simplifies our analysis, and existing upper bounds in the literature are derived in the same setting. Our contributions are as follows:

• For SGD with random reshuffling, we prove a lower bound of . We conjecture that it is tight, as it seems to combine the “best” behaviors of existing upper bounds: It behaves as for a small constant number

of passes (which is optimal as discussed above), interpolating to

when is large enough, and contains a term decaying cubically with . Moreover, the lower bound holds under more general conditions than the upper bounds: For example, it holds for any , and even if the function under consideration is quadratic and on .

• For SGD with a single shuffling, we prove a lower bound of . Although we are not aware of an upper bound to compare to, this lower bound already proves an inherent performance gap compared to random reshuffling: Indeed, in the latter case there is an upper bound of , which is smaller than the lower bound for single shuffling when is sufficiently large. This implies that the added computational effort of repeatedly reshuffling the functions can provably pay off in terms of the convergence rate.

• For cyclic SGD, we provide an lower bound. We note that a similar bound (at least asymptotically and for a certain ) is already implied by [5, Theorem 3.4]. Our contribution here is to present a more explicit and non-asymptotic lower bound which holds for any and .

## 2 Preliminaries

We let bold-face letters denote vectors. A twice-differentiable function on is -strongly convex, if its Hessian satisfies for all . is quadratic if it is of the form for some matrix , vector and scalar .

We consider finite-sum optimization problems as in Eq. (1), and our lower bound constructions satisfy the following rather specific conditions for some given positive parameters (recall that for lower bounds, making the assumptions more stringent actually strengthens the result):

###### Assumption 1.

is a quadratic finite-sum function of the form for some , which is -strongly convex. Each is convex and quadratic, has -Lipschitz gradients, and moreover, is -Lipschitz for any such that where . Also, the algorithm is initialized at some for which .

Before continuing, we make a few remarks about the setting and our results:

###### Remark 1 (Unconstrained Optimization).

For simplicity, in this paper we consider unconstrained SGD, where the iterates are not explicitly constrained to lie in some subset of the domain. However, we note that existing upper bounds for SGD on strongly convex functions often assume an explicit projection on such a subset, in order to ensure that the gradients remain bounded. That being said, it is not difficult to verify that all our constructions – which have a very simple structure – are such that the iterates remain in a region with bounded gradients (with probability

, at least for reasonably small step sizes), in which case projections will not significantly affect the results.

###### Remark 2 (Distance from Optimum).

In Assumption 1, we fix the initial distance from the optimum to be at most , rather than keeping it as a variable parameter. Besides simplifying the constructions, we note that existing SGD upper bounds for strongly convex functions often do not explicitly depend on the initial distance (both for with-replacement SGD and with random reshuffling, see for example [10, 12, 8]). Thus, it makes sense to study lower bounds in which the initial distance is fixed to be some constant.

###### Remark 3 (Applicability of the Lower Bounds).

We emphasize that in our lower bounds, we focus on (a) SGD with constant step size, and (b) the expected performance of the iterate after exactly epochs. Thus, they do not formally cover step sizes which change across iterations, the performance of other iterates, or the performance of some average of the iterates. However, it is not clear that these are truly necessary to achieve optimal error bounds (indeed, many existing analyses do not require them), and we conjecture that our lower bounds cannot be substantially improved even with non-constant step sizes and iterate averaging schemes.

## 3 SGD with Random Reshuffling

We begin by discussing SGD with random reshuffling, where at the beginning of every epoch we choose a new random order for processing the individual functions (Algorithm 1). Our main result is the following:

###### Theorem 1.

For any , and positive such that , there exists a function on and an initialization point satisfying Assumption 1, such that for any step size ,

 E[F(xk)−infxF(x)] ≥ c⋅min{λ , G2λ(1(nk)2+1nk3)} ,

where is a universal constant.

For large, this lower bound is

 Ω(G2λ(1(nk)2+1nk3)) .

It is useful to compare this bound to the existing optimal bound for SGD with replacement, which is

 Θ(G2λnk)

(see for example ). First, we note that the factor is the same in both of them. The dependence on though is different: For or constant , our lower bound is , similar to the with-replacement case, but as increases, it decreases cubically (rather than linearly) with . This indicates that even for small , random reshuffling is superior to with-replacement sampling, which agrees with empirical observations. For very large (

), a phase transition occurs and the bound becomes

– that is, scaling down quadratically with the total number of individual stochastic iterations. That being said, it should be emphasized that is often an unrealistic regime, especially in large-scale problems where is a huge number.

The proof of Thm. 1 appears in Sec. 6.3. It is based on a set of very simple constructions, where , and the individual functions are all of the form for appropriate . This allows us to write down the iterates at the end of each epoch in closed form. The analysis then carefully tracks the decay of after each epoch, showing that it cannot decay to too rapidly, hence implying a lower bound on after epochs. The main challenge is that unlike SGD with replacement, here the stochastic iterations in each epoch are not independent, so computing these expectations is not easy. To make it tractable, we identify two distinct sources contributing to the error in each epoch: A “bias” term, which captures the fact that the stochastic gradients at each epoch are statistically correlated, hence for a given iterate during the algorithm’s run,

(unlike the with-replacement case where equality holds), and a “variance” term, which captures the inherent noise in the stochastic sampling process. For different parameter regimes, we use different constructions and focus on either the bias or the variance component (which when studied in isolation are more tractable), and then combine the various bounds into the final lower bound appearing in Thm.

1.

We finish with the following remark about a possible extension of the lower bound:

###### Remark 4 (Convex Functions).

For convex functions which are not necessarily strongly convex (namely, any is possible), Thm. 1 seems to suggest a lower bound (in terms of ) of

 Ω⎛⎝G√1(nk)2+1nk3⎞⎠ = Ω(G(1√nk3+1nk)) ,

since in this scenario we can set arbitrarily small, and in particular as so as to maximize the lower bound in Thm. 1. In contrast,  show a upper bound in this setting for SGD with random reshuffling, and a similar upper bound hold for SGD with replacement. A similar argument can also be applied to the other lower bounds in our paper, extending them from the strongly convex to the convex case. However, we emphasize that some caution is needed, since our lower bounds do not quantify a dependence on the radius of the domain, which is usually explicit in bounds for this setting. We leave a formal lower bound result in the general convex case for future work.

## 4 SGD with a Single Shuffling

We now turn to the case of SGD where a single random order over the individual functions is chosen at the beginning, and the algorithm then cycles over that order (Algorithm 2). Our main result is the following:

###### Theorem 2.

For any , and positive such that , there exists a function on and an initialization point satisfying Assumption 1, such that for any step size ,

 E[F(xk)−infxF(x)] ≥ c⋅min{λ , G2λnk2} ,

where is a universal constant.

The proof appears in Subsection 6.2. In the single shuffling case, we are not aware of a good upper bound to compare to (except the bound for cyclic SGD below, which trivially applies also the SGD with single shuffling). However, the lower bound already implies an interesting separation between single shuffling and random reshuffling: In the former case, is the best we can hope to achieve, whereas in the latter case, we have seen upper bounds which are strictly better when is sufficiently large (i.e., ). To the best of our knowledge, this is the first formal separation between these two shuffling schemes for SGD: It implies that the added computational effort of repeatedly reshuffling the functions can provably pay off in terms of the convergence rate. It would be quite interesting to understand whether this separation might also occur for smaller values of as well, which is definitely true if our lower bound for random reshuffling is tight. It would also be interesting to derive a good upper bound for SGD with single shuffling, which is a common heuristic.

## 5 Cyclic SGD

Finally, we turn to discuss cyclic SGD, where the individual functions are cycled over in a fixed deterministic order. We note that for this algorithm, an lower bound was already proven in , but in an asymptotic form, and only for . Our contribution here is to provide an explicit bound which holds for any :

###### Theorem 3.

For any , and positive such that , there exists a function on and an initialization point satisfying Assumption 1, such that if we run cyclic GD for epochs with any step size , then

 F(xk)−infxF(x) ≥ c⋅min{λ , G2λk2}

where is a universal constant.

The proof (which follows a broadly similar strategy to Thm. 1) appears in Sec. 6.3. Comparing this theorem with our other lower bounds and the associated upper bounds, it is clear that there is a high price to pay (in a worst-case sense) for using a fixed, non-random order, as the bound does not improve at all with more individual functions . Indeed, recalling that the bound for with-replacement SGD is , it follows that cyclic SGD can beat with-replacement SGD only when , or . For large-scale problems where is big, this often an unrealistically large number of passes.

## 6 Proofs

### 6.1 Proof of Thm. 1

For simplicity, we will prove the theorem assuming the number of components in our function is an even number. This is without loss of generality, since if

is odd, let

be the function achieving the lower bound using an even number of components, and define where . has the same Lipschitz parameter as , and a strong convexity parameter smaller than that of by a factor which is always in . Moreover, it is easy to see that for a fixed step size, the distribution of the iterates after epochs is the same over and , since SGD does not move on any iteration where is chosen. Therefore, the lower bound on translates to a lower bound on up to a small factor which can be absorbed into the numerical constants. Thus, in what follows, we will assume that is even and that , whereas in the theorem statement we make the slightly stronger assumption so that the reduction described above will be valid.

The proof of the theorem is based on the following three propositions, each using a somewhat different construction and analysis:

###### Proposition 1.

For any even and any positive such that , there exists a function on satisfying Assumption 1, such that for any step size ,

 E[F(xk)−infxF(x)] ≥ c⋅min{λ , G2λnk3}

where is a universal constant.

###### Proposition 2.

Suppose that and that is even. For any positive such that , there exists a function on satisfying Assumption 1, such that for any step size ,

 E[F(xk)−infxF(x)] ≥ c⋅G2λ(nk)2

where is a numerical constant.

###### Proposition 3.

Suppose and that is even. For any positive such that , there exists a function on satisfying Assumption 1, such that for any step size ,

 E[F(xk)−infxF(x)] ≥ c⋅min{λ , G2λ(nk)2}

where is a numerical constant.

The proof of each proposition appears below, but let us first show how combining these implies our theorem. We consider two cases:

• If , then , so by Proposition 1,

 E[F(xk)−infxF(x)] ≥ c⋅min{λ , G2λnk3} ≥ c⋅min{λ , G22λ(1(nk)2+1nk3)} .
• If (which implies since is even), we have , and by combining Proposition 2 and Proposition 3 (which together cover any positive step size),

 E[F(xk)−infxF(x)] ≥ c⋅min{λ , G2λ(nk)2} ≥ c⋅min{λ , G22λ(1(nk)2+1nk3)}

Thus, in any case we get , from which the result follows.

#### 6.1.1 Proof of Proposition 1

We will need the following key technical lemma, whose proof (which is rather long and technical) appears in Appendix A:

###### Lemma 1.

Let (for even ) be a random permutation of (where both and appear exactly times). Then there is a numerical constant , such that for any ,

 E⎡⎣(n−1∑i=0σi(1−α)i)2⎤⎦ ≥ c⋅min{1+1α , n3α2}

Let be fixed (assuming and is even). We will use the following function:

 F(x)=1nn∑i=1fi(x) = λ2x2 ,

where , and

 fi(x)={λ2x2+G2xi≤n2λ2x2−G2xi>n2 . (2)

Also, we assume that the algorithm is initialized at . On this function, we have that during any single epoch, we perform iterations of the form

 xnew=(1−ηλ)xold+ηG2σi,

where are a random permutation of ’s and ’s. Repeatedly applying this inequality, we get that after iterations, the relationship between the first and last iterates in the epoch satisfy

 xt+1 = (1−ηλ)nxt+ηG2n−1∑i=0σi(1−ηλ)n−i−1 = (1−ηλ)nxt+ηG2n−1∑i=0σi(1−ηλ)i . (3)

(in the last equality, we used the fact that are exchangeable). Using this and the fact that , we get that

 E[x2t+1] = (1−ηλ)2nE[x2t]+(ηG2)2⋅βn,η,λ , (4)

where

 βn,η,λ := E⎡⎣(n−1∑i=0σi(1−λη)n−i−1)2⎤⎦ = E⎡⎣(n−1∑i=0σi(1−λη)i)2⎤⎦. (5)

Note that if , then by Lemma 1, for some positive constant , and we get that

 E[x2t+1] ≥ (ηG2)2⋅c ≥ (G2λ)2⋅c

for all , and therefore , so the proposition we wish to prove holds. Thus, we will assume from now on that .

With this assumption, repeatedly applying Eq. (4) and recalling that , we have

 E[x2k] ≥ (1−ηλ)2nk+(ηG2)2⋅βn,η,λk−1∑t=0(1−ηλ)2nt = (1−ηλ)2nk+(ηG2)2⋅βn,η,λ⋅1−(1−ηλ)2nk1−(1−ηλ)2n . (6)

We now consider a few cases (recalling that the case was already treated earlier):

• If , then we have

 E[x2k] ≥ (1−ηλ)2nk ≥ (1−12nk)2nk≥14

for all .

• If then by Bernoulli’s inequality, we have , and therefore, by Eq. (6)

 E[x2k] ≥ η2G2βn,η,λ(1−(1−1/2nk)2nk)4(1−(1−2nηλ)) ≥ ηG2βn,η,λ(1−exp(−1))8λn .

Plugging in Lemma 1 and simplifying a bit, this is at least

 cηG2λn⋅min{1ηλ,n3(ηλ)2} = cηG2λn⋅n3(ηλ)2 = cη3λn2G2

for some numerical constant . Using the assumption that (which implies ), this is at least

 c8⋅G2λ2nk3 .
• If , then is at least some numerical constant , so Eq. (6) implies

 E[x2k] ≥ c(ηG2)2⋅βn,η,λ .

By Lemma 1, this is at least

 c′(ηG2)2⋅min{1+1ηλ , n3(ηλ)2} = c′(ηG2)2(1+1ηλ) ≥ c′ηG24λ

Since , this is at least

 c′G28λ2n ≥ c′G28λ2nk3 .

Combining all the cases, we get overall that

 E[x2k] ≥ c⋅min{1,G2λ2nk3}

for some numerical constant . Noting that and combining with the above, the result follows.

#### 6.1.2 Proof of Proposition 2

We use the same construction as in the proof of Proposition 1, where , and leading to Eq. (6), namely

 E[x2k] ≥ (1−ηλ)2nk+(ηG2)2⋅βn,η,λ⋅1−(1−ηλ)2nk1−(1−ηλ)2n , (7)

where , are a random permutation of ’s and ’s.

As in the proof of Proposition 1, we consider several regimes of . In the same manner as in that proof, it is easy to verify that when or , then is at least a positive constant (hence ) , and when , for a numerical constant (hence ). In both these cases, the statement in our proposition follows, so it is enough to consider the regime .

In this regime, by Bernoulli’s inequality, we have , so we can lower bound Eq. (7) by

 (ηG2)2⋅βn,η,λ1−(1−ηλ)2nk2nηλ = ηG2βn,η,λ(1−(1−ηλ)2nk)8λn .

Since we assume , it follows that for some positive . Plugging this and the bound for from Lemma 1, the displayed equation above is at least

 cηG28λn⋅min{1ηλ,n3(ηλ)2} = cηG28λn⋅n3(ηλ)2 = c8G2λη3n2 .

Since we assume , this is at least

 c′⋅G2λ2n4

for some numerical . Since we assume that , this is at least . Noting that and combining with the above, the result follows.

#### 6.1.3 Proof of Proposition 3

To simplify some of the notation, we will prove the result for a function which is -strongly convex (rather than -strongly convex), assuming , and notice that this only affects the universal constant in the bound. Specifically, we use the following function:

 F(x)=1nn∑i=1fi(x) = λ4x2 ,

where , and

 fi(x)={λ2x2+G2xi≤n2−G2xi>n2 .

Also, we assume that the algorithm is initialized at . On this function, we have that during any single epoch, we perform iterations of the form

 xnew=(1−ηλσi)xold+ηG2(1−2σi),

where are a random permutation of ’s and ’s. Repeatedly applying this equation, we get that after iterations, the relationship between the iterates and is

 xt+1 = xt⋅n−1∏i=0(1−ηλσi)+ηG2n−1∑i=0(1−2σi)n−1∏j=i+1(1−ηλσj) (8)

As a result, and using the fact that are independent of and in , we have

 E[x2t+1] ≥ E[x2t⋅n−1∏i=0(1−ηλσi)2]+ηG⋅E[xt(n−1∏i=0(1−ηλσi))(n−1∑i=0(1−2σi)n−1∏j=i+1(1−ηλσj))] ≥ (1−ηλ)2n⋅E[x2t]+ηG⋅E[xt]⋅E[(n−1∏i=0(1−ηλσi))(n−1∑i=0(1−2σi)n−1∏j=i+1(1−ηλσj))] (9)

We now wish to use Lemma 6 from Appendix B, in order to replace the products in the expression above by sums. To that end, and in order to simplify the notation, define

 A:=n−1∏i=0(1−ηλσi)  ,  Bi:=n−1∏j=i+1(1−ηλσj)  ,  ~A:=1−ηλn∑i=1σi=1−ηλn2  ,  ~Bi:=1−ηλn∑j=i+1σi , (10)

and note that by Lemma 6,

 An−1∑i=0(1−2σi)Bi ≤ ⎛⎝~A±2(ηλn−1∑i=0σi)2⎞⎠⎛⎝n−1∑i=0(1−2σi)~Bi±2n−1∑i=0(ηλn−1∑j=i+1σj)2⎞⎠ , (11)

where is taken to be either plus or minus depending on the sign of and , to make the inequality valid (we note that eventually we will show that these terms are relatively negligible). Opening the product, and using the deterministic upper bounds

 |~A|≤1  ,  (ηλn−1∑i=0σi)2≤(ηλn)2 (12)

and

 ∣∣ ∣∣n−1∑i=0(1−2σi)~Bi∣∣ ∣∣≤n  ,  n−1∑i=0(ηλn−1∑j=i+1σj)2 ≤ n(ηλn)2 ≤1104n , (13)

(which follow from the assumption that ), we can upper bound Eq. (11) by

 ~An−1∑i=0(1−2σi)~Bi+2(ηλn)2⋅(n+2100n)+n(ηλn)2 (∗)≤ ~An−1∑i=0(1−2σi)~Bi+301100(ηλ)2n3 ,

where in we used the fact that and therefore . Substituting back the definitions of and plugging back into Eq. (11), we get that

 E [(n−1∏i=0(1−ηλσi))(n−1∑i=0(1−2σi)n−1∏j=i+1(1−ηλσj))] (∗)≤ ηλn(−(1−ηλn2)n+14(n−1)+301100ηλn2) ,

where is by Lemma 4. Using the assumptions that (hence ) and , this is at most for a numerical constant . Summarizing this part of the proof, we have shown that

 E[(n−1∏i=0(1−ηλσi))(n−1∑i=0