Proceedings of Machine Learning Research vol 336:1–24, 2026
39th Annual Conference on Learning Theory
Random Reshuffling Dominates Stochastic Gradient Descent Zijian Liu
ZL 3067@ STERN . NYU . EDU
Stern School of Business, New York University
arXiv:2606.32005v1 [math.OC] 30 Jun 2026
Editors: Steve Hanneke and Tor Lattimore
Abstract Stochastic Gradient Descent (SGD) is one of the most classical optimization algorithms with favorable theoretical guarantees, yet the practical implementation of SGD differs subtly from its well-known form and is often referred to as Shuffling Stochastic Gradient Descent (Shuffling SGD). A particularly popular strategy in Shuffling SGD is Random Reshuffling (RR), which has achieved great empirical success across numerous experiments. Despite its strong performance, RR has long been considered a heuristic due to a lack of theoretical support. Over the last decade, people have finally established provable convergence rates for RR, thus justifying its observed superiority. However, for smooth convex optimization, two clouds over the convergence theory of RR remain to this day. More precisely, according to the current theory, Shuffling SGD under RR converges only when the stepsize is smaller than a threshold proportional to 1/n, where n is the number of summands in the objective (or the number of data points). Consequently, the optimally tuned theoretical rate of Shuffling SGD under RR is strictly worse than that of SGD when the number of epochs is smaller than another threshold proportional to n. These two restrictions heavily limit the applicability of existing theories and leave a critical mismatch with practice. In this work, for the first time, we prove that RR dominates SGD in smooth convex optimization under any reasonable stepsize after any finite number of epochs, thereby addressing a longstanding open question. Keywords: Convex Optimization, Stochastic Optimization, Random Reshuffling
1. Introduction One of the fundamental tasks in machine learning is to optimize functions in a finite-sum form, P i.e., f (x) ≜ n1 ni=1 fi (x). Among different optimization algorithms, Stochastic Gradient Descent (SGD), proposed in the seminal work of Robbins and Monro (1951), is arguably one of the most classical methods. Due to its easy implementation and computational efficiency, SGD is particularly popular when n is large, the standard case nowadays. More importantly, the convergence guarantees of SGD have been extensively studied, yielding provable rates in various settings (Polyak, 1987; Bottou et al., 2018; Lan, 2020), thereby providing a theoretical backbone for SGD. However, compared with the standard form of SGD analyzed in theory, which uniformly samples a function to perform a gradient descent step at each iteration, the practical implementation differs subtly and is often referred to as Shuffling Stochastic Gradient Descent (Shuffling SGD). In Shuffling SGD, the optimization procedure is divided into K epochs, and within each epoch, the order in which functions are processed is determined by a permutation π of {1, . . . , n}. A widely implemented strategy for generating π is Random Reshuffling (RR), which, in each epoch, independently and uniformly draws a new permutation from all possible ones. Although Shuffling SGD under RR has achieved great empirical success across numerous experiments, it has long been considered a heuristic due to a lack of theoretical support. Over 0. In this work, we say that one optimization algorithm dominates another if the order of its convergence rate is no worse than that of the latter and is strictly better in certain regimes.
© 2026 Z. Liu.
L IU
the last decade, beginning with the pioneering work of Gürbüzbalaban et al. (2021), people have finally established provable convergence rates for RR, thereby justifying its observed superiority over standard SGD. In particular, for smooth convex optimization (i.e., each fi is convex and L-smooth1 ), RR with a D2 constant stepsize2 η is known to converge in expectation at the rate ηnK +η 2 nLσ⋆2 (e.g., Mishchenko 1 et al. (2020); Nguyen et al. (2021)), provided that the stepsize satisfies η ≲ nL , where D denotes 2 the distance between the initial point and the optimal solution, and σ⋆ is the gradient variance at the minimizer. In comparison, SGD under the same setting guarantees the in-expectation convergence D2 + ησ⋆2 (Garrigos and Gower, 2023) but only requires η ≲ L1 . Clearly, RR converges faster rate ηnK 1 than SGD in the regime η ≲ nL , which has been recognized as theoretical evidence demonstrating the strong performance of RR. Despite the progress discussed above, some important issues remain unaddressed. The most critical longstanding open question is that people still do not understand what happens to RR when 1 η falls into the regime η ≳ nL . This point is critical because, given that n is typically large in 1 modern tasks, the threshold nL can be extremely small or even vanish, whereas stepsizes used in practice are usually at a constant level, thereby leaving a significant gap between theory and practice. More crucially, even if one temporarily assumes that the existing rate for RR mentioned earlier could be extended to allow η ≲ L1 (though no such theory has been established), it would still fail to explain the advantage of RR over SGD, as the term η 2 nLσ⋆2 for RR is worse than the 1 term ησ⋆2 for SGD when η ≳ nL . This hints that an analysis different from existing ones may be needed. Another issue implied by the above discussion is that the optimally tuned rate for RR induced 2 Lσ⋆2 D4 1 σ⋆ D LD2 √ 3 by the existing bound is only LD K + ( nK 2 ) , which is better than the best tuned rate nK + nK for SGD only when K is larger than a threshold proportional to n. Moreover, in the case of σ⋆ = 0 2 (i.e., all fi ’s share a common optimal solution), the rate of RR reduces to only LD K , which is even 2 1 worse than the LD nK rate of SGD by a factor of n . The above restrictions on the stepsize η or the number of epochs K heavily limit the applicability of existing theories and cannot fully explain the favorable performance of Shuffling SGD under RR compared with standard SGD. Therefore, we are naturally led to the following question: In smooth convex optimization, does RR dominate SGD without these two restrictions? 1.1. Our Contributions This work provides an affirmative answer to the above question. 2
D • Concretely, we show that Shuffling SGD under RR provably converges at a rate of ηnK + 1 2 min {1, ηnL} ησ⋆ for any stepsize satisfying η ≲ L (see Theorem 1 for the formal version with nonuniform smoothness parameters and dynamic stepsizes that depend on the epoch number). 1 We highlight that this rate is not only the first provable result for RR that allows η ≳ nL , but also D2 1 2 provably dominates the ηnK + ησ⋆ bound of SGD under any stepsize η ≲ L . It is noteworthy that this rate is not merely a simple extension of the previously best bound for RR, since the latter 1 is slower than SGD when η ≳ nL , as discussed before.
1. For simplicity, we adopt a uniform smoothness parameter in the discussion, as in most of the existing literature. 2. We also use a constant stepsize in the discussion for convenience.
2
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
2
• Consequently, the optimally tuned rate for Shuffling SGD under RR is improved to LD nK + 2
4
1
D ⋆D , ( Lσ min{ √σ⋆nK ) 3 } (see Corollary 1 for the formal version with nonuniform smoothness parameters), nK 2 2
σ⋆ D √ which dominates the best tuned bound LD nK + nK of SGD for any finite K. Moreover, the rate 2
2
LD 1 reduces to LD nK when σ⋆ = 0, improving upon the best known result K by a factor of n .
In summary, for the first time, we prove that RR dominates standard SGD in smooth convex optimization under any reasonable stepsize after any finite number of epochs, resolving a longstanding open question. 1.2. Related Work We provide a brief overview of Shuffling SGD under RR and defer further details to Appendix A. Over the past decades, the effectiveness of RR has been reported in many works (e.g., Bottou (2009, 2012); Bengio (2012)). However, the theoretical understanding of it has long lagged behind. The first breakthrough is by Gürbüzbalaban et al. (2021), which provides the first theoretical evidence that RR can beat SGD in smooth strongly convex optimization under certain additional assumptions. Since then, RR has been extensively studied (e.g., Nagaraj et al. (2019); Haochen and Sra (2019)). D2 To date, the best known rate in smooth convex optimization is ηnK + η 2 nLσ⋆2 under the condition 1 η ≲ nL (Mishchenko et al., 2020; Nguyen et al., 2021). The only existing lower bound in smooth 2
4
1
⋆D convex optimization is ( Lσ ) 3 , due to Cha et al. (2023), which holds for constant stepsizes nK 2 2 2 1 η ≲ nL and a large number of epochs at least satisfying K ≳ nLσ2D . ⋆
2. Preliminary Notation. N denotes the set of natural numbers (excluding 0). Given p n ∈ N, we write [n] ≜ {1, . . . , n}. ⟨·, ·⟩ is the standard Euclidean inner product, and ∥·∥ ≜ ⟨·, ·⟩ is the ℓ2 norm. Given a real-valued differentiable function h : Rd → R, ∇h(x) denotes the gradient at x ∈ Rd . The Bregman divergence induced by h is defined as Bh (x, y) ≜ h(x) − h(y) − ⟨∇h(y), x − y⟩, which is nonnegative if h is additionally convex. Objective.
We study the following finite-sum optimization problem in this work n
1X inf f (x) ≜ fi (x), n x∈Rd i=1
where n ∈ N and each fi : Rd → R is differentiable. Remark 1 To ease notation, hereinafter we use B ≜ Bf and Bi ≜ Bfi to denote the Bregman divergences induced by f and fi , respectively. Assumptions.
Our analysis relies on the following three assumptions.
Assumption 1 (Minimizer) ∃x⋆ ∈ Rd such that f⋆ ≜ f (x⋆ ) = inf x∈Rd f (x) ∈ R. Assumption 2 (Convexity) Each fi is convex. 3
L IU
Assumption 3 (Smoothness) Each fi is Li -smooth, i.e., ∃Li > 0 such that ∥∇fi (x) − ∇fi (y)∥ ≤ Li ∥x − y∥ , ∀x, y ∈ Rd . All three of the above assumptions are standard and commonly adopted in the literature (Polyak, 1987; Nesterov et al., 2018; Bottou et al., 2018; Lan, 2020). Notably, we do not impose any assumptions on the difference between the individual gradient ∇fi and the full gradient ∇f , such as the popular finite variance condition. Next, we introduce three more notations to simplify the expressions in the subsequent sections. σ⋆2 denotes the variance of the gradient at the minimizer x⋆ , and L̄ (resp. L̂) represents the average (resp. maximum) smoothness parameter, i.e., n
σ⋆2 ≜
P 1X ∥∇fi (x⋆ )∥2 , L̄ ≜ n1 ni=1 Li , L̂ ≜ max Li . n i∈[n] i=1
We note that the quantity σ⋆2 is widely used in prior works on shuffling gradient methods (e.g., Ying et al. (2019); Mishchenko et al. (2020); Nguyen et al. (2021)) and remains invariant even when f has multiple minimizers (see Lemma 4.17 of Garrigos and Gower (2023)). In particular, σ⋆2 = 0 corresponds to the case in which all fi share a common optimal solution. To finish this section, we state a classical result in convex optimization, known as the cocoercivity property of smooth convex functions, which serves as a key tool in our analysis. As for its proof, see, for example, Theorem 2.15 of Nesterov et al. (2018). Lemma 1 (Co-coercivity) Let h : Rd → R be a differentiable convex function that is also Lsmooth, then we have, for any x, y ∈ Rd , ∥∇h(x) − ∇h(y)∥2 ≤ 2LBh (x, y) and ∥∇h(x) − ∇h(y)∥2 ≤ L ⟨∇h(x) − ∇h(y), x − y⟩ .
3. Shuffling Stochastic Gradient Descent Algorithm 1 Shuffling Stochastic Gradient Descent (Shuffling SGD) Input: initial point x11 ∈ Rd , stepsize ηk > 0 for k = 1 to K do Generate a permutation πk of [n] for i = 1 to n do xi+1 = xik − ηk ∇fπi (xik ) k k end for x1k+1 = xn+1 k end for The method studied in this work, Shuffling Stochastic Gradient Descent (Shuffling SGD), is given in Algorithm 1. Compared with the standard SGD algorithm, which uniformly samples a function to process at each step, Shuffling SGD determines the order in which functions are passed in each epoch based on a permutation. In particular, three strategies for generating permutations are popular in practice, as illustrated in the following examples, where Sn denotes the symmetric group of [n]. Example 1 (Random Reshuffling (RR)) Each πk is drawn independently and uniformly from Sn . 4
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
Example 2 (Single Shuffling (SS)) Each πk = π, a permutation drawn uniformly from Sn . Example 3 (Incremental Gradient (IG)) Each πk = π, a deterministic permutation from Sn . 3.1. New Rate for RR We are now ready to provide the main result, Theorem 1, a new convergence rate for RR. Theorem 1 Under Assumptions 1, 2, and 3, suppose RR is employed with ηk ≤ 1 , ∀k ∈ [K], let 6L̂
x̄K ≜
K X n K X X ηk i xk where HK ≜ ηk , nHK k=1 i=1
(1)
k=1
then Shuffling SGD (Algorithm 1) guarantees that 2
6 x11 − x⋆ + 51 min E [f (x̄K ) − f⋆ ] ≤ P η n K k k=1
(P
PK K 2 3 k=1 ηk k=1 ηk nL̄ , PK P K k=1 ηk k=1 ηk
) σ⋆2 .
Remark 2 We make no effort to optimize the constants in the bounds obtained in this work. Proof The proof is deferred to Subsection 4.4. To the best of our knowledge, in smooth convex optimization, Theorem 1 shows the first convergence rate for RR under any reasonable stepsize, i.e., ηk ≲ 1/L̂, thereby improving existing 2 results in different aspects. Previously, the best known convergence rate for RR is PD + n K k=1 ηk PK p 2 3 k=1 ηk nL̄σ⋆ PK (where D ≜ x11 − x⋆ ) under the condition ηk ≲ 1/(n L̄L̂) (Liu and Zhou, 2024; k=1 ηk
Cai et al., 2024). However, this requirement on ηk is highly unrealistic, as n is typically pretty large in modern machine learning tasks, meaning that the stepsize has to be extremely small or even vanishing. This contradicts the constant-level stepsizes commonly used in practice and leads to a gap between theory and practice. In comparison, Theorem 1 not only allows the stepsize to lie in the constant regime ηk ≲ 1/L̂ but also shows a fundamental improvement rather than a mere extension of the previously best known rate, since the latter becomes slower than SGD once ηk ≳ 1/(nL̄) (the rate of SGD is stated below), while our Theorem 1 never does. P In addition, compared with the
n
2 PD K
k=1 ηk
K
+
ηk2 σ⋆2 k=1 ηk
k=1 P K
rate of standard SGD (Garrigos and
Gower, 2023), our Theorem 1 is never worse for any reasonable stepsize (i.e., ηk ≲ 1/L̂), and is strictly better when ηk ≲ 1/(nL̄). This feature has an important implication: if both methods employ their own optimally tuned stepsizes, Shuffling SGD under RR provably achieves a better 2 σ⋆ D √ upper bound that dominates the L̂D nK + nK rate of SGD after any finite number of epochs, as evidenced by Corollary 1 below. Corollary 1 Under the same setting as in Theorem 1, with the optimally tuned constant stepsize ηk = η⋆ , ∀k ∈ [K], where η⋆ ≤ 1 , Shuffling SGD (Algorithm 1) guarantees that 6L̂ !1 2 σ x1 − x 2 x1 − x 4 3 L̂ x11 − x⋆ L̄σ ⋆ ⋆ ⋆ ⋆ 1 √1 E [f (x̄K ) − f⋆ ] ≲ + min , . 2 nK nK nK 5
L IU
Proof The proof is deferred to Subsection 4.4. √ 2 4 1 2 ⋆D In contrast, the known optimally tuned rate for RR is only L̄KL̂D + ( L̄σ ) 3 . To better nK 2 understand the differences among these optimally tuned rates, let OrderSGD (K) and OrderOld RR (K) denote the dominant terms in the existing optimally tuned rates for SGD and Shuffling SGD under RR, respectively, i.e., ( ) (p 1 ) L̂D2 σ⋆ D L̄L̂D2 L̄σ⋆2 D4 3 Old OrderSGD (K) ≜ max . ,√ and OrderRR (K) ≜ max , nK K nK 2 nK Similarly, OrderNew RR (K) denotes the dominant term in the rate obtained in Corollary 1, i.e., ( ( 2 4 31 )) L̂D2 σ⋆ D L̄σ⋆ D New OrderRR (K) ≜ max . , min √ , nK nK 2 nK In the nondegenerate case σ⋆ ̸= 0 (which implies that n ≥ 2), the following comparisons hold: (a) New nL̄2 D2 , (K) = Order (K) ≤ OrderOld Order SGD RR RR (K), K ≤ σ⋆2 (b) 2 Old nL̄2 D2 OrderNew < K ≤ nL̄σL̂D , 2 RR (K) < OrderSGD (K) ≤ OrderRR (K), σ⋆2 ⋆ (c) 1/2 3/2 2 Old nL̄L̂D2 OrderNew < K ≤ nL̄ σL̂2 D , RR (K) ≤ OrderRR (K) < OrderSGD (K), σ⋆2 ⋆ OrderNew (K) = OrderOld (K) < Order nL̄1/2 L̂3/2 D2 (K), K > , 2 SGD RR RR σ ⋆
2
2
where (a) becomes an equality if and only if K = nL̄σ2D and L̄ = L̂, (b) becomes an equality if ⋆
1/2
2
3/2
2
, and (c) becomes an equality if and only if K = nL̄ σL̂2 D . As one can and only if K = nL̄σL̂D 2 ⋆ ⋆ see, unlike the existing optimally tuned bound for RR in the literature, which can be slower than 2 standard SGD when K ≤ nL̄σL̂D , the new result in Corollary 1 is never worse than the optimally 2 ⋆
2
2
tuned rate of standard SGD and strictly beats it once K > nL̄σ2D , thereby improving the threshold ⋆
2
2
2
from nL̄σL̂D to nL̄σ2D . 2 ⋆ ⋆ Moreover, in the special case where all fi ’s share a common minimizer, or equivalently when 2 rate as SGD. However, the prior bound for RR can only σ⋆ = 0, Corollary 1 achieves the same L̂D nK √ p 2 be reduced to a slower rate in the order of L̄KL̂D due to the stepsize restriction ηk ≲ 1/(n L̄L̂), as discussed before. Formally, once n ≥ 2, we have Old OrderNew RR (K) = OrderSGD (K) < OrderRR (K), ∀K ∈ N.
In summary, Theorem 1 and Corollary 1 together imply that RR dominates SGD in smooth convex optimization under any reasonable stepsize after any finite number of epochs.
4. Theoretical Analysis In this section, we lay the groundwork for proving Theorem 1 and complete its proof at the end. The section is organized into four parts. First, we provide the two most important lemmas in Subsection 6
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
4.1. Next, in Subsection 4.2, we establish an upper bound in Theorem 2, which indicates that Shuffling SGD under RR never converges more slowly than SGD under any reasonable stepsize. Then, Theorem 3 in Subsection 4.3 presents an alternative convergence rate for Shuffling SGD under RR, which demonstrates that Algorithm 1 under RR provably converges faster than SGD when the stepsize is sufficiently small. Finally, in Subsection 4.4, we conclude Theorem 1 and then use it to prove Corollary 1. 4.1. Two Core Lemmas This subsection contains two core lemmas, both of which are critical to our analysis. Before presenting the lemmas, we introduce two notions. Given a permutation π of [n] and two indices i, j ∈ [n] satisfying j ≤ i, we define the following new permutation π(i, j) ≜ π 1 , . . . , π j−1 , π i , π j+1 , . . . , π i−1 , π j , π i+1 , . . . , π n .
(2)
In words, π(i, j) is the permutation generated by exchanging the elements πi and πj in π. Equipped with the notion of π(i, j), we introduce the following virtual sequence, for any given k ∈ [K], 1 1 l l xl+1 k (i, j) ≜ xk (i, j) − ηk ∇fπ l (i,j) (xk (i, j)), ∀l ∈ [n] , where xk (i, j) ≜ xk . k
(3)
This means that the sequence xlk (i, j), ∀l ∈ [n + 1] denotes the trajectory of the k-th epoch starting from x1k , produced by Algorithm 1, but under the permutation πk (i, j). This virtual iterate can be viewed as a coupled sequence of the real output and, to the best of our knowledge, was first introduced by Sherman et al. (2021). It plays a fundamental role in our proof, as will become clear. With these two concepts in hand, we proceed to state the two core lemmas. The first is Lemma 2, which is based on Lemma 2 of Sherman et al. (2021). For completeness, we reproduce the proof of Lemma 2 in Appendix C. P Lemma 2 Given an arbitrary finite-sum function ℓ(x) = n1 ni=1 ℓi (x), suppose RR is employed, then for any k ∈ [K] and i ∈ [n], Shuffling SGD (Algorithm 1) guarantees that i 1X h i h E ℓ(xik ) − ℓπi (xik ) = E ℓπi (xik (i, j)) − ℓπi (xik ) , k k k n j<i
where xk (i, j) is defined in (3). In the analysis of shufflingigradient methods, a well-known major challenge is to properly bound the h i term E f (xk ) − fπi (xik ) , unlike in SGD, which no longer equals 0 due to the nature of shufflingk based algorithms. Lemma 2 provides a possible approach by relating the term we want to control (in a slightly more general form, applicable to any finite-sum function ℓ) to another quantity involving the virtual sequence introduced earlier in (3). For the convenience of the discussion, temporarily assume ℓi = fi in Lemma 2. Then, under the smoothness assumption, one would expect the difference between ℓπi (xik (i, j)) and ℓπi (xik ) to k k be small whenever xik (i, j) and xik are close. This observation naturally leads us to the other core Lemma 3 stated below. 7
L IU
Lemma 3 Under Assumptions 2 and 3, suppose ηk ≤ 2 , ∀k ∈ [K], then for any k ∈ [K], i ∈ [n], L̂ and j ∈ [i − 1], Shuffling SGD (Algorithm 1) guarantees that xik (i, j) − xik ≤ ηk ∇fπi (xjk ) − ∇fπj (xjk ) , k
k
where xk (i, j) is defined in (3). Lemma 3 quantifies how close the virtual iterate and the true trajectory can be, under the widely required condition of ηk ≤ 2/L̂ in smooth optimization. Although the inequality does not directly offer a bound on the distance between xik (i, j) and xik that depends only on deterministic terms (e.g., the stepsize ηk ), it is sufficient for our proof when combined with a careful analysis. The proof of Lemma 3 is given in Appendix C and relies on the co-coercivity property (i.e., Lemma 1), which is closely related to the nonexpansiveness of the update rule in gradient methods (Polyak, 1987; Nesterov et al., 2018). 4.2. Bound I: Never Worse than SGD under Reasonable Stepsize In this subsection, we give the first bound for Shuffling SGD under RR, stated in Theorem 2 below. Theorem 2 Under Assumptions 1, 2, and 3, suppose RR is employed with ηk ≤ 1 , ∀k ∈ [K], 6L̂ then Shuffling SGD (Algorithm 1) guarantees that P 2 2 2 x11 − x⋆ 51 K k=1 ηk σ⋆ + E [f (x̄K ) − f⋆ ] ≤ , P PK 2n K k=1 ηk k=1 ηk where x̄K is defined in (1). Discussion on Theorem 2. To the best of our knowledge, Theorem 2 offers the first theoretical evidence that Shuffling SGD under RR shares surprising similarities with SGD, as reflected in the two aspects elaborated below. First, Theorem 2 states that, similar to SGD, Shuffling SGD under RR does converge under any reasonable stepsize (i.e., ηk ≲ 1/L̂). In contrast, as far as we know, all prior works that provide provable rates of multi-epoch RR for smooth convex optimization require the stepsize ηk to be smaller than a threshold proportional to 1/n, with only one exception (Nagaraj et al., 2019), which, however, assumes each fi to be additionally Lipschitz, thereby limiting the applicability of their theory and even excluding common quadratic optimization problems over Rd . Second, we highlight that Theorem 2 gives the same convergence upper bound (up to constant factors) as SGD (Garrigos and Gower, 2023) in smooth convex optimization, while allowing a stepsize that depends on the epoch number. This result thus fills a gap in the literature. Putting these together, Theorem 2 indicates that, for smooth convex optimization, RR under any reasonable stepsize never converges more slowly than SGD. Analysis. In the following, we present the analysis for Theorem 2 and finally prove it. The core idea underlying the proof is, as one might expect, to analyze Shuffling SGD in a manner analogous to SGD. In other words, we aim to quantify the progress made by Algorithm 1 at each iteration. Although this perspective is natural, it has been less explored in prior studies. The main reason is that, as discussed earlier, the permutation in Shuffling SGD causes the most important property of 8
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
SGD, unbiasedness, to no longer hold. To overcome this barrier, we develop a novel analysis that avoids any additional assumptions, such as Lipschitz continuity considered in Nagaraj et al. (2019). We start with the following Lemma 4, a standard step in characterizing the per-iterate progress of Shuffling SGD (or SGD). The proof of Lemma 4 follows directly from expanding both sides. To make the work self-contained, we include it in Appendix D. Lemma 4 Under Assumption 1, for any k ∈ [K] and i ∈ [n], Shuffling SGD (Algorithm 1) guarantees that fπi (xik ) − fπi (x⋆ ) = k k
xik − x⋆
2
2
− xi+1 − x⋆ k 2ηk
+
2 ηk ∇fπi (xik ) − Bπi (x⋆ , xik ). k k 2
By the co-coercivity property (i.e., Lemma 1), the term η2k ∥∇fπi (xik )∥2 −Bπi (x⋆ , xik ) can be easily k
k
upper bounded by ηk ∥∇fπi (x⋆ )∥2 (ignoring constant factors) once ηk ≲ 1/L̂, which further yields k a desired residual term ηk σ⋆2 after taking expectations. Therefore, the only difficulty is to relate fπi (xik ) to f (xik ), which is, again, the main challenge in the analysis of Shuffling SGD. k To address the issue mentioned, the prior work of Nagaraj et al. (2019) applies an argument based on Wasserstein distance, which additionally requires the Lipschitz continuity of each fi . In comparison, we tackle this problem by establishing the following new inequality in Lemma 5. Lemma 5 Under Assumptions 2 and 3, suppose RR is employed with ηk ≤ 1 , ∀k ∈ [K], then for 6L̂ any k ∈ [K] and i ∈ [n], Shuffling SGD (Algorithm 1) guarantees that E
f (xik )
i h 2 2 4ηk X j i i + ≤ E fπi (xk ) + ηk E ∇fπi (xk ) . E ∇fπj (xk ) k k k 3n j<i
Remark 3 Lemma 5 is stronger than the existing bound of Nagaraj et al. (2019) derived via the Wasserstein distance, since imposing the additional condition ∥∇fi (x)∥ ≤ G, as in Nagaraj et al. (2019), recovers their Lemma 4. Lemma 5 provides a novel inequality that measures the difference between E f (xik ) and E[fπi (xik )] k by the second moment of the stochastic gradients up to time i. Note that the second term on the R.H.S. can be absorbed by the R.H.S. of the inequality in Lemma 4. For the remaining term, the coefficient ηk /n is P key to the final proof, which ensures that the accumulated error in one epoch is controlled by E[ηk ni=1 ∥∇fπi (xik )∥2 ]. k The proof of Lemma 5 builds on the two core results, Lemmas 2 and 3, presented before. To save space, we defer it to Appendix D. Final proof. With Lemmas 4 and 5 stated above, we are ready to prove Theorem 2. Proof of Theorem 2 We sum the inequality in Lemma 4 from i = 1 to n and use x1k+1 = xn+1 to k obtain n X i=1
fπi (xik )−fπi (x⋆ ) = k k
x1k − x⋆
2
− x1k+1 − x⋆ 2ηk 9
2
+
n
n
i=1
i=1
2 X ηk X ∇fπi (xik ) − Bπi (x⋆ , xik ). k k 2
L IU
h i Take expectations on both sides and note that E fπi (x⋆ ) = f⋆ to yield k
n X
h h
i
E fπi (xik ) − f⋆ =
2
x1k − x⋆
E
i
−E
h
2
x1k+1 − x⋆
i
2ηk
k
i=1
+E
" n X ηk
∇fπi (xik )
2
i=1
2
# − Bπi (x⋆ , xik ) . k
k
Next, we invoke Lemma 5 and sum it up from i = 1 to n to have # " n n n h i X 2 X 7ηk X i i i E f (xk ) ≤ . E fπi (xk ) + E ∇fπi (xk ) k k 3 i=1
i=1
(4)
(5)
i=1
Combine (4) and (5) to obtain n X
h E E f (xik ) − f⋆ ≤
x1k − x⋆
2
i
−E
h
x1k+1 − x⋆
2
i
2ηk
i=1
+E
" n X 17ηk i=1
6
∇fπi (xik ) k
2
# − Bπi (x⋆ , xik ) . k
(6)
One more step, we observe that 2 2 18 ∇fπi (xik ) − ∇fπi (x⋆ ) + 18 ∇fπi (x⋆ ) k k k 17 (a) 36 2 ≤ Lπi Bπi (x⋆ , xik ) + 18 ∇fπi (x⋆ ) k 17" k k " # # n n X X (b) 2 17ηk ⇒E ≤E Bπi (x⋆ , xik ) + 51ηk nσ⋆2 , ∇fπi (xik ) k k 6
∇fπi (xik ) k
i=1
2
≤
(7)
i=1
where (a) is by Lemma 1 and (b) is due to ηk ≤ 1 and E 6L̂
∇fπi (x⋆ ) k
2
= σ⋆2 , ∀k ∈ [K] , i ∈ [n].
Finally,Pwe plug (7) back into (6), multiply both sides by ηk , sum over k = 1 to K, divide both sides by n K k=1 ηk , apply the convexity of f , and use the definition of x̄K in (1) to conclude. 4.3. Bound II: Always Better than SGD under Small Stepsize This subsection presents the other rate of Shuffling SGD under RR, as shown in Theorem 3 below. Theorem 3 Under Assumptions 1, 2, and 3, suppose RR is employed with ηk ≤ 1 , ∀k ∈ [K], 2L̂ then Shuffling SGD (Algorithm 1) guarantees that P 2 3 2 6 x11 − x⋆ 8 K k=1 ηk nL̄σ⋆ E [f (x̄K ) − f⋆ ] ≤ + , PK P K n k=1 ηk k=1 ηk where x̄K is defined in (1). 10
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
Discussion on Theorem 3. Readers familiar with the literature on shuffling gradient methods may readily figure out that the rate given in Theorem 3 perfectly matches the known bound for Shuffling SGD under RR in smooth convex optimization (e.g., Mishchenko et al. (2020); Nguyen et al. (2021)). However, we emphasize a key difference here, that is, the stepsize in our Theorem 3 is allowed to satisfy ηk ≲ 1/L̂, in contrast to all existing results that require ηk to be at most inversely proportional to n. More importantly, in the setting of nonuniform Li considered in this work, the largest threshold on the p stepsize in the literature that guarantees a rate similar to Theorem 3 is in the order of 1/(n L̄L̂) (Liu and Zhou, 2024; Cai et al., 2024). But as indicated by our Theorem 3, the superiority of RR over SGD already exists once ηk ≲ 1/(nL̄). Especially, this improvement can be significant when a dominant smoothness parameter exists, leading to L̂ ≈ nL̄. Therefore, Theorem 3 is the first result to extend the known bound in smooth convex optimization to any reasonable stepsize while preserving the favorable property of Algorithm 1, i.e., Shuffling SGD under RR provably converges faster than SGD when the stepsize is sufficiently small. Analysis. The roadmap for establishing Theorem 3 differs wildly from that used before to prove Theorem 2. This time, our proof strategy is to check how close Shuffling SGD can be to Gradient Descent. More concretely, we will view each epoch of Algorithm 1 (containing n iterations) as a single step and analyze the progress made by it at once. This kind of approach has appeared in different previous works (e.g., Mishchenko et al. (2020); Nguyen et al. (2021)) and always yields a convergence rate in the order of 1/(n1/3 K 2/3 ), faster than SGD when K is large. However, all works that follow the idea described above share the same issue, that is, they require the stepsize to be in the order of 1/n, which is, however, not enough for our purpose. In the following, we develop a new analysis to bypass this critical obstacle. We now formally begin the analysis by introducing another virtual sequence, defined as follows, for any k ∈ [K], 1 ≜ ykn+1 , where y11 ≜ x⋆ . yki+1 ≜ yki − ηk ∇fπi (x⋆ ), ∀i ∈ [n] , yk+1 k
(8)
Under the above definition, and noting that ∇f (x⋆ ) = 0, one can find ykn+1 = yk1 − ηk
n X
∇fπi (x⋆ ) = yk1 − ηk n∇f (x⋆ ) = yk1 , ∀k ∈ [K] . k
i=1
Combine the above line and y11 = x⋆ as defined in (8) to have ykn+1 = yk1 = x⋆ , ∀k ∈ [K] .
(9)
The above virtual sequence is inspired by the work of Mishchenko et al. (2020), which, as far as we know, was the first to propose a similar term under the constant stepsize. Here, we slightly extend their idea to accommodate the case where the stepsize can depend on the current epoch number. Remark 4 We note that Mishchenko et al. (2020) introduced the virtual sequence to handle the case of individual strong convexity, i.e., each fi is required to be strongly convex. However, in our setting, only individual convexity is assumed. This means that their proof cannot be applied. As such, our analysis substantially departs from the existing approach. 11
L IU
Equipped with the new virtual sequence introduced above, we first present the following Lemma 6. Lemma 6 Under Assumptions 1, 2, and 3, suppose ηk ≤ 1 , ∀k ∈ [K], then for any k ∈ [K], 2L̂ Shuffling SGD (Algorithm 1) guarantees that ηk
n X
Bπi (xik , x⋆ ) ≤ x1k − x⋆ k
2
− x1k+1 − x⋆
i=1
2
+ 2ηk
n X
Bπi (yki , x⋆ ), k
i=1
where yki is defined in (8). As discussed earlier, we intentionally treat each epoch of Shuffling SGD as a single step. Hence, compared with Lemma 4 used to prove the first bound, Lemma 6 is in a different flavor, which shows the progress made by Algorithm 1 over an entire epoch. two tasks naturally P arise. The first is to lower bound the term PnBased oni the form of Lemma 6, P n n i,x ) = i B (x , x ) on the L.H.S. by B(x i ⋆ ⋆ i=1 πk i=1 i=1 f (xk ) − f⋆ . The second is to upper k k Pn bound the residual term i=1 Bπi (yki , x⋆ ) on the R.H.S. k We address the first task by establishing Lemma h7, a novel inequality upper i P h the following i j 1 i i 2 bounding each E B(xk , x⋆ ) in terms of E Bπi (xk , x⋆ ) , n j<i E Bπj (xk , x⋆ ) , and σ⋆ . k
k
Lemma 7 Under Assumptions 1, 2, and 3, suppose RR is employed with ηk ≤ 1 , ∀k ∈ [K], then 2L̂ for any k ∈ [K] and i ∈ [n], Shuffling SGD (Algorithm 1) guarantees that i i 4X h h E Bπj (xjk , x⋆ ) + 12ηk2 (i − 1)L̄σ⋆2 . E B(xik , x⋆ ) ≤ 2E Bπi (xik , x⋆ ) + k k n j<i
To the best of our knowledge, we are the first to obtain an inequality in such a hform. Intuitively, i Lemma 7 says that each E B(xik , x⋆ ) differs from its stochastic counterpart E Bπi (xik , x⋆ ) by hk i at most a multiplicative constant, together with an average of the preceding terms E Bπj (xjk , x⋆ ) k satisfying j < i (not an exact average due to the coefficient 4/n), plus an additional term involving σ⋆2h. Summing the inequality in Lemma 7 from i = 1 to n yields a meaningful lower bound on i Pn Pn 2 i i E i=1 B(xk , x⋆ ) and σ⋆ . i=1 Bπki (xk , x⋆ ) based on E The proof of Lemma 7 is rather technical, so we skip the discussion here. For details, we kindly refer the interested reader to Appendix E. P Lastly, we need to bound the residual term ni=1 Bπi (yki , x⋆ ) in Lemma 6. The inequality k obtained is given in Lemma 8 below. Lemma 8 Under Assumptions 1 and 3, suppose RR is employed, then for any k ∈ [K], Shuffling SGD (Algorithm 1) guarantees that " n # X η 2 n2 L̄σ⋆2 E Bπi (yki , x⋆ ) ≤ k , k 6 i=1
where yki is defined in (8). Lemma 8 can be derived in a relatively easy way, as it can be deduced from existing works. In particular, thanks to Liu and Zhou (2024), we obtain a bound that depends only on the average smoothness parameter L̄ rather than the maximum smoothness parameter L̂. 12
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
Final proof. Armed with Lemmas 6, 7, and 8 above, we are finally able to prove Theorem 3. Proof of Theorem 3 First, we invoke Lemma 6 to have ηk
n X
n X 2 1 xk+1 − x⋆ + 2ηk Bπi (yki , x⋆ ), k i=1
2 x1k − x⋆ −
Bπi (xik , x⋆ ) ≤ k
i=1
where yki is defined in (8). Take expectations on both sides and apply Lemma 8 to yield ηk E
" n X
# Bπi (xik , x⋆ ) ≤ E
h
k
x1k − x⋆
2
i
−E
h
x1k+1 − x⋆
2
i
ηk3 n2 L̄σ⋆2 . 3
+
i=1
(10)
Next, we multiply both sides of the inequality in Lemma 7 by ηk and sum it up from i = 1 to n to obtain " n # " n # X X ηk E (11) B(xik , x⋆ ) ≤ 6ηk E Bπi (xik , x⋆ ) + 6ηk3 n2 L̄σ⋆2 . k
i=1
i=1
Combine (10) and (11) to have ηk E
" n X
# B(xik , x⋆ )
≤ 6E
h
2
x1k − x⋆
i
− 6E
h
x1k+1 − x⋆
2
i
+ 8ηk3 n2 L̄σ⋆2
i=1
⇒E
"K n XX
# ηk B(xik , x⋆ ) ≤ 6 x11 − x⋆
2
+
k=1 i=1
K X
8ηk3 n2 L̄σ⋆2 .
k=1
Finally, we observe that B(xik , x⋆ ) = f (xik ) − f⋆ , divide both sides by n convexity of f , and use the definition of x̄K in (1) to conclude.
PK
k=1 ηk , apply the
4.4. Proofs of Theorem 1 and Corollary 1 With the previous preparation, Theorem 1 follows immediately. Proof of Theorem 1 Combine Theorems 2 and 3 to conclude. We next derive Corollary 1 directly from Theorem 1. Proof of Corollary 1 With a constant stepsize ηk = η ≤ 1 , ∀k ∈ [K], Theorem 1 reduces to 6L̂
6 x11 − x⋆ E [f (x̄K ) − f⋆ ] ≤ ηnK
2
+ 51 min η, η 2 nL̄ σ⋆2 .
Optimizing the R.H.S. of the above inequality over 0 < η ≤ 1 yields the desired result. 6L̂
13
L IU
5. Conclusion and Future Work In this work, we prove that Shuffling SGD under RR dominates SGD in smooth convex optimization under any reasonable stepsize after any finite number of epochs. Our main Theorem 1 follows from combining two novel convergence results, whose analysis may each be of independent interest. Our work suggests several new directions for future research. From an upper-bound perspective, our current proof is split into two distinct parts. It is therefore worthwhile to investigate whether Theorem 1 can be obtained via a unified analysis. From a lower-bound perspective, the only existing hardness result for RR in smooth convex optimization by Cha et al. (2023) is established for the stepsize ηk at most inversely proportional to n and for the number of epochs K at least proportional to n. As such, providing a complete characterization of the lower bound for RR under any reasonable stepsize and any number of epochs remains an important task for future work.
Acknowledgments The author thanks the anonymous reviewers for their valuable feedback.
References Kwangjun Ahn, Chulhee Yun, and Suvrit Sra. Sgd with shuffling: optimal rates without component convexity and large epoch requirements. In H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin, editors, Advances in Neural Information Processing Systems, volume 33, pages 17526–17535. Curran Associates, Inc., 2020. URL https://proceedings.neurips.cc/paper_files/paper/2020/ file/cb8acb1dc9821bf74e6ca9068032d623-Paper.pdf. Yoshua Bengio. Practical Recommendations for Gradient-Based Training of Deep Architectures, pages 437–478. Springer Berlin Heidelberg, Berlin, Heidelberg, 2012. ISBN 978-3-64235289-8. doi: 10.1007/978-3-642-35289-8_26. URL https://doi.org/10.1007/ 978-3-642-35289-8_26. Léon Bottou. Curiously fast convergence of some stochastic gradient descent algorithms. In Proceedings of the symposium on learning and data science, Paris, volume 8, pages 2624–2633. Citeseer, 2009. Léon Bottou. Stochastic Gradient Descent Tricks, pages 421–436. Springer Berlin Heidelberg, Berlin, Heidelberg, 2012. ISBN 978-3-642-35289-8. doi: 10.1007/978-3-642-35289-8_25. URL https://doi.org/10.1007/978-3-642-35289-8_25. Léon Bottou, Frank E. Curtis, and Jorge Nocedal. Optimization methods for large-scale machine learning. SIAM Review, 60(2):223–311, 2018. doi: 10.1137/16M1080173. URL https:// doi.org/10.1137/16M1080173. Xufeng Cai and Jelena Diakonikolas. Last iterate convergence of incremental methods as a model of forgetting. In Y. Yue, A. Garg, N. Peng, F. Sha, and R. Yu, editors, International Conference on Learning Representations, volume 2025, pages 102613–102647, 2025. URL https://proceedings.iclr.cc/paper_files/paper/2025/file/ fea9f93f4cec99f65a8b4d575fc353a8-Paper-Conference.pdf. 14
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
Xufeng Cai, Cheuk Yin Lin, and Jelena Diakonikolas. Tighter convergence bounds for shuffled sgd via primal-dual perspective. In A. Globerson, L. Mackey, D. Belgrave, A. Fan, U. Paquet, J. Tomczak, and C. Zhang, editors, Advances in Neural Information Processing Systems, volume 37, pages 72475–72524. Curran Associates, Inc., 2024. doi: 10.52202/ 079017-2310. URL https://proceedings.neurips.cc/paper_files/paper/ 2024/file/84d395725a9b40cb4a49d84478ac24c7-Paper-Conference.pdf. Jaeyoung Cha, Jaewook Lee, and Chulhee Yun. Tighter lower bounds for shuffling SGD: Random permutations and beyond. In Andreas Krause, Emma Brunskill, Kyunghyun Cho, Barbara Engelhardt, Sivan Sabato, and Jonathan Scarlett, editors, Proceedings of the 40th International Conference on Machine Learning, volume 202 of Proceedings of Machine Learning Research, pages 3855–3912. PMLR, 23–29 Jul 2023. URL https://proceedings.mlr.press/ v202/cha23a.html. Guillaume Garrigos and Robert M Gower. Handbook of convergence theorems for (stochastic) gradient methods. arXiv preprint arXiv:2301.11235, 2023. Mert Gürbüzbalaban, Asu Ozdaglar, and Pablo A Parrilo. Why random reshuffling beats stochastic gradient descent. Mathematical Programming, 186:49–84, 2021. Jeff Haochen and Suvrit Sra. Random shuffling beats SGD after finite epochs. In Kamalika Chaudhuri and Ruslan Salakhutdinov, editors, Proceedings of the 36th International Conference on Machine Learning, volume 97 of Proceedings of Machine Learning Research, pages 2624–2633. PMLR, 09–15 Jun 2019. URL https://proceedings.mlr.press/v97/ haochen19a.html. Vladimir Kibardin. Decomposition into functions in the minimization problem. Automation and Remote Control, 1979, 01 1979. Tomer Koren, Roi Livni, Yishay Mansour, and Uri Sherman. Benign underfitting of stochastic gradient descent. In S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh, editors, Advances in Neural Information Processing Systems, volume 35, pages 19605–19617. Curran Associates, Inc., 2022. URL https://proceedings.neurips.cc/paper_files/paper/2022/file/ 7bc4f74e35bcfe8cfe43b0a860786d6a-Paper-Conference.pdf. Guanghui Lan. First-order and stochastic optimization methods for machine learning. Springer, 2020. Zijian Liu and Zhengyuan Zhou. On the last-iterate convergence of shuffling gradient methods. In Ruslan Salakhutdinov, Zico Kolter, Katherine Heller, Adrian Weller, Nuria Oliver, Jonathan Scarlett, and Felix Berkenkamp, editors, Proceedings of the 41st International Conference on Machine Learning, volume 235 of Proceedings of Machine Learning Research, pages 32471– 32508. PMLR, 21–27 Jul 2024. URL https://proceedings.mlr.press/v235/ liu24cg.html. Zijian Liu and Zhengyuan Zhou. Improved last-iterate convergence of shuffling gradient methods for nonsmooth convex optimization. In Aarti Singh, Maryam Fazel, Daniel Hsu, Simon LacosteJulien, Felix Berkenkamp, Tegan Maharaj, Kiri Wagstaff, and Jerry Zhu, editors, Proceedings 15
L IU
of the 42nd International Conference on Machine Learning, volume 267 of Proceedings of Machine Learning Research, pages 40152–40193. PMLR, 13–19 Jul 2025. URL https: //proceedings.mlr.press/v267/liu25ct.html. Konstantin Mishchenko, Ahmed Khaled, and Peter Richtarik. Random reshuffling: Simple analysis with vast improvements. In H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin, editors, Advances in Neural Information Processing Systems, volume 33, pages 17309–17320. Curran Associates, Inc., 2020. URL https://proceedings.neurips.cc/paper_files/paper/2020/file/ c8cc6e90ccbff44c9cee23611711cdc4-Paper.pdf. Dheeraj Nagaraj, Prateek Jain, and Praneeth Netrapalli. SGD without replacement: Sharper rates for general smooth convex functions. In Kamalika Chaudhuri and Ruslan Salakhutdinov, editors, Proceedings of the 36th International Conference on Machine Learning, volume 97 of Proceedings of Machine Learning Research, pages 4703–4711. PMLR, 09–15 Jun 2019. URL https://proceedings.mlr.press/v97/nagaraj19a.html. Angelia Nedic and Dimitri P. Bertsekas. Incremental subgradient methods for nondifferentiable optimization. SIAM Journal on Optimization, 12(1):109–138, 2001. doi: 10.1137/ S1052623499362111. URL https://doi.org/10.1137/S1052623499362111. Yurii Nesterov et al. Lectures on convex optimization, volume 137. Springer, 2018. Lam M. Nguyen, Quoc Tran-Dinh, Dzung T. Phan, Phuong Ha Nguyen, and Marten van Dijk. A unified convergence analysis for shuffling-type gradient methods. Journal of Machine Learning Research, 22(207):1–44, 2021. URL http://jmlr.org/papers/v22/20-1238.html. Boris T. Polyak. Introduction to optimization. New York, Optimization Software, 1987. Shashank Rajput, Anant Gupta, and Dimitris Papailiopoulos. Closing the convergence gap of SGD without replacement. In Hal Daumé III and Aarti Singh, editors, Proceedings of the 37th International Conference on Machine Learning, volume 119 of Proceedings of Machine Learning Research, pages 7964–7973. PMLR, 13–18 Jul 2020. URL https://proceedings.mlr. press/v119/rajput20a.html. Shashank Rajput, Kangwook Lee, and Dimitris Papailiopoulos. Permutation-based SGD: Is random optimal? In International Conference on Learning Representations, 2022. URL https:// openreview.net/forum?id=YiBa9HKTyXE. Herbert Robbins and Sutton Monro. A Stochastic Approximation Method. The Annals of Mathematical Statistics, 22(3):400 – 407, 1951. doi: 10.1214/aoms/1177729586. URL https: //doi.org/10.1214/aoms/1177729586. Itay Safran and Ohad Shamir. How good is sgd with random shuffling? In Jacob Abernethy and Shivani Agarwal, editors, Proceedings of Thirty Third Conference on Learning Theory, volume 125 of Proceedings of Machine Learning Research, pages 3250–3284. PMLR, 09–12 Jul 2020. URL https://proceedings.mlr.press/v125/safran20a.html. 16
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
Itay Safran and Ohad Shamir. Random shuffling beats sgd only after many epochs on ill-conditioned problems. In M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang, and J. Wortman Vaughan, editors, Advances in Neural Information Processing Systems, volume 34, pages 15151–15161. Curran Associates, Inc., 2021. URL https://proceedings.neurips.cc/paper_files/paper/2021/file/ 803ef56843860e4a48fc4cdb3065e8ce-Paper.pdf. Ohad Shamir. Without-replacement sampling for stochastic gradient methods. In D. Lee, M. Sugiyama, U. Luxburg, I. Guyon, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 29. Curran Associates, Inc., 2016. URL https://proceedings.neurips.cc/paper_files/paper/2016/file/ c74d97b01eae257e44aa9d5bade97baf-Paper.pdf. Uri Sherman, Tomer Koren, and Yishay Mansour. Optimal rates for random order online optimization. In M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang, and J. Wortman Vaughan, editors, Advances in Neural Information Processing Systems, volume 34, pages 2097–2108. Curran Associates, Inc., 2021. URL https://proceedings.neurips.cc/paper_ files/paper/2021/file/107030ca685076c0ed5e054e2c3ed940-Paper. pdf. Bicheng Ying, Kun Yuan, Stefan Vlaski, and Ali H. Sayed. Stochastic learning under random reshuffling with constant step-sizes. IEEE Transactions on Signal Processing, 67(2):474–489, 2019. doi: 10.1109/TSP.2018.2878551.
17
L IU
Appendix A. Additional Related Work This section provides additional discussion of the related work. We mainly focus on smooth convex optimization under RR. As for nonsmooth convex optimization under RR/SS/IG, the reader could refer to Kibardin (1979); Nedic and Bertsekas (2001); Koren et al. (2022); Liu and Zhou (2025). See also Shamir (2016) for RR under structured problems. The first breakthrough on RR is by Gürbüzbalaban et al. (2021) for smooth strongly convex optimization, which, however, requires each fi to be quadratic or to have a Lipschitz Hessian. Since then, extensive studies have emerged. Among them, a series of works continues to study the convergence behavior of RR in smooth strongly convex optimization for quadratic objectives or under other additional conditions (Ying et al., 2019; Haochen and Sra, 2019; Safran and Shamir, 2020; Rajput et al., 2020; Ahn et al., 2020; Safran and Shamir, 2021; Rajput et al., 2022). To the best of our knowledge, the first work that drops the strong convexity assumption is D2 Nagaraj et al. (2019), which provides a convergence rate of ηnK + ηG2 under the requirement η ≤ L2 for L-smooth G-Lipschitz convex fi , where we remind the reader that D denotes the distance between the initial point and the optimal solution, η represents the stepsize, and K is the number of epochs. However, this rate cannot reflect any advantage of RR over standard SGD. Subsequently, two works (Mishchenko et al., 2020; Nguyen et al., 2021) further remove the extra D2 1 Lipschitz assumption and establish the bound ηnK + η 2 nLσ⋆2 under the condition η ≲ nL , where 2 σ⋆ is the gradient variance at the optimal solution. This rate is faster than standard SGD under its 1 required regime and remains the best bound so far. In fact, it is unimprovable for small η ≲ nL 2
2
2
4
1
⋆D and large K ≳ nLσ2D due to the lower bound of ( Lσ ) 3 by Cha et al. (2023). Recently, Liu and nK 2 ⋆ Zhou (2024) and Cai and Diakonikolas (2025) extend the above rate from the average iterate to the last iterate (up to additional polylogarithmic factors).
Appendix B. Summary of Notation For readability, we recall and summarize the notation used in the paper. P • f = n1 ni=1 fi is the objective, where each fi : Rd → R is differentiable and convex. P • x⋆ denotes the minimizer of f . f⋆ = f (x⋆ ) is the optimal function value. σ⋆2 = n1 ni=1 ∥∇fi (x⋆ )∥2 is the variance of the gradient at the optimal solution. P • Li > 0 is the smoothness parameter of fi . L̄ = n1 ni=1 Li is the average smoothness parameter. L̂ = maxi∈[n] Li is the maximum smoothness parameter. • Bh (x, y) = h(x) − h(y) − ⟨∇h(y), x − y⟩ denotes the Bregman divergence induced by a realvalued differentiable function h (not necessarily convex). We write B (resp. Bi ) as shorthand for Bf (resp. Bfi ).
Appendix C. Missing Proofs of Core Lemmas This section contains the missing proofs of the two core lemmas presented in Subsection 4.1. Before providing the proofs, we recall two notions introduced in (2) and (3), respectively. Given a permutation π of [n] and two indices i, j ∈ [n] satisfying j ≤ i, π(i, j) is the permutation 18
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
generated by exchanging the elements πi and πj in π, i.e., π(i, j) = π 1 , . . . , π j−1 , π i , π j+1 , . . . , π i−1 , π j , π i+1 , . . . , π n .
(12)
For any given k ∈ [K], xlk (i, j), ∀l ∈ [n + 1] denotes the trajectory of the k-th epoch starting from x1k , produced by Shuffling SGD, but under the permutation πk (i, j), i.e., 1 1 l l xl+1 k (i, j) = xk (i, j) − ηk ∇fπ l (i,j) (xk (i, j)), ∀l ∈ [n] , where xk (i, j) = xk . k
(13)
C.1. Proof of Lemma 2 Proof Given k ∈ [K] and i ∈ [n], we first have the decomposition
ℓ(xik ) − ℓπi (xik ) = k
1 n
X j<i
ℓπj (xik ) − ℓπi (xik ) +
X
k
k
ℓπj (xik ) − ℓπi (xik ) . k
j≥i
k
i i h h Note that E ℓπj (xik ) = E ℓπi (xik ) holds for any j ∈ {i, . . . , n} under RR, since ℓπj (xik ) and k
k
k
ℓπi (xik ) are equal in distribution conditioning on πk1 to πki−1 and π1 to πk−1 . Therefore, we obtain k
i 1X h h i E ℓ(xik ) − ℓπi (xik ) = E ℓπj (xik ) − ℓπi (xik ) . k k k n
(14)
j<i
For any fixed j ∈ [i], it is known that πk equals πk (i, j) in distribution (e.g., Lemma C.3 of Liu and Zhou (2025)), which implies that (πk , xk ) also equals (πk (i, j), xk (i, j)) in distribution, since the trajectory of the k-th epoch generated by Shuffling SGD is deterministically determined by the permutation and stepsize. This implies that i i h i h h i i i E ℓπj (xk ) = E ℓπj (i,j) (xk (i, j)) = E ℓπi (xk (i, j)) , k
k
k
where the last step is due to πkj (i, j) = πki from its definition (2). Hence, we finally obtain h i 1X h i E ℓπi (xik (i, j)) − ℓπi (xik ) . E ℓ(xik ) − ℓπi (xik ) = k k k n j<i
C.2. Proof of Lemma 3 Proof By the definition of π(i, j) in (2), we know πkl (i, j) = πkl , ∀l ∈ / {i, j} ⇒ ∇fπl (i,j) = ∇fπl , ∀l ∈ / {i, j} . k
19
k
(15)
L IU
Therefore, given l ∈ / {i, j}, by the definition of xk (i, j) (see (3)) and the update rule of Shuffling SGD, we have 2 2 l+1 l l l l xl+1 (i, j) − x = x (i, j) − x − η ∇f (x (i, j)) − ∇f (x ) l l k k k k πk (i,j) k πk k k 2 (15) = xlk (i, j) − xlk − ηk ∇fπl (xlk (i, j)) − ∇fπl (xlk ) k
k
2
2
= xlk (i, j) − xlk + ηk2 ∇fπl (xlk (i, j)) − ∇fπl (xlk ) k k D E l l l − 2ηk ∇fπl (xk (i, j)) − ∇fπl (xk ), xk (i, j) − xlk k k ! (a) 2 2 2ηk 2 l l ≤ xk (i, j) − xk + ηk − ∇fπl (xlk (i, j)) − ∇fπl (xlk ) k k Lπl k
(b)
2
≤ xlk (i, j) − xlk
,
(16)
where (a) is due to Lemma 1 and (b) holds by ηk ≤ 2 ⇒ ηk2 − L2ηkl ≤ 0. L̂
π
Apply (16) from l = j + 1 to l = i − 1 to obtain
k
2
2
j+1 ≤ xj+1 k (i, j) − xk 2 = xjk (i, j) − xjk − ηk ∇fπj (i,j) (xjk (i, j)) − ∇fπj (xjk )
xik (i, j) − xik
k
k
(c) 2 =ηk
∇fπj (i,j) (xjk ) − ∇fπj (xjk )
2 (2)
= ηk2 ∇fπi (xjk ) − ∇fπj (xjk ) k
k
k
2
,
k
(3)
where (c) holds by xjk (i, j) = xjk , since xlk (i, j) = xlk , ∀l ∈ [j] implied by (16) and x1k (i, j) = x1k together.
Appendix D. Missing Proofs of Lemmas for Bound I In this section, we provide the missing proofs of the lemmas presented in Subsection 4.2, which were used to prove the first convergence rate in Theorem 2. D.1. Proof of Lemma 4 Proof Given k ∈ [K] and i ∈ [n], by the definition of Bπi , k
D
E
fπi (xik ) − fπi (x⋆ ) = ∇fπi (xik ), xik − x⋆ − Bπi (x⋆ , xik ) k
k
k
=
xik − x⋆
k
2
− xi+1 − x⋆ k 2ηk
2
+
2 ηk ∇fπi (xik ) − Bπi (x⋆ , xik ), k k 2
where the second step holds by the update rule of Algorithm 1.
20
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
D.2. Proof of Lemma 5 Proof We apply Lemma 2 with ℓ = f and ℓi = fi to have h i 1X h i E f (xik ) − fπi (xik ) = E fπi (xik (i, j)) − fπi (xik ) , k k k n
(17)
j<i
where xk (i, j) is defined in (3). Next, by the Lπi -smoothness of fπi (Assumption 3), we know k
k
fπi (xik (i, j)) − fπi (xik ) k k D E Lπi 2 ≤ ∇fπi (xik ), xik (i, j) − xik + k xik (i, j) − xik k 2 (a) ηk2 Lπi 2 k ≤ ηk ∇fπi (xik ) ∇fπi (xjk ) − ∇fπj (xjk ) + ∇fπi (xjk ) − ∇fπj (xjk ) k k k k k 2 2 (b) ηk + 2ηk Lπi 2 2 k ≤ηk ∇fπi (xik ) + ∇fπi (xjk ) − ∇fπj (xjk ) k k k 4 ηk + 2ηk2 Lπi 2 2 2 j j i k ≤ηk ∇fπi (xk ) + ∇fπi (xk ) + ∇fπj (xk ) k k k 2 (c) 2 2 2 2 ≤ηk ∇fπi (xik ) + ηk , ∇fπi (xjk ) + ∇fπj (xjk ) k k k 3
(18)
where (a) is by Cauchy-Schwarz inequality and Lemma 3, (b) is due to AM-GM inequality, and (c) holds by ηk ≤ 1 . Finally, we plug (18) back into (17) to obtain 6L̂ i η X h 2 2 2 2 2 k j j i i i E f (xk ) − fπi (xk ) ≤ E ∇fπi (xk ) + ∇fπi (xk ) + ∇fπj (xk ) k k k k n 3 3 j<i 2 2 4 (d) ηk X = E ∇fπi (xik ) + ∇fπj (xjk ) k k n 3 j<i 2 2 4ηk X j i ≤ ηk E ∇fπi (xk ) + E ∇fπj (xk ) , k k 3n j<i
where (d) is due to E
∇fπi (xjk ) k
2
=E
∇fπj (xjk ) k
2
when j < i.
Appendix E. Missing Proofs of Lemmas for Bound II In this section, we provide the missing proofs of the lemmas presented in Subsection 4.3, which were used to prove the second convergence rate in Theorem 3. Before presenting the proofs, we recall the notion introduced in (8). For any k ∈ [K], the virtual sequence yki , ∀i ∈ [n + 1] follows the equations, 1 yki+1 = yki − ηk ∇fπi (x⋆ ), ∀i ∈ [n] , yk+1 = ykn+1 , where y11 = x⋆ . k
(19)
As shown in (9), the virtual sequence satisfies that ykn+1 = yk1 = x⋆ , ∀k ∈ [K] . 21
(20)
L IU
E.1. Proof of Lemma 6 Proof Given k ∈ [K] and i ∈ [n], by the update rule of Algorithm 1 and the definition of yki+1 in (8), we have 2 = xik − yki − ηk ∇fπi (xik ) − ∇fπi (x⋆ ) k k D E 2 2 xik − yki + ηk2 ∇fπi (xik ) − ∇fπi (x⋆ ) + 2ηk ∇fπi (x⋆ ) − ∇fπi (xik ), xik − yki k k k k 2 2 xik − yki + ηk2 ∇fπi (xik ) − ∇fπi (x⋆ ) + 2ηk Bπi (yki , x⋆ ) − Bπi (yki , xik ) − Bπi (xik , x⋆ ) k k k k k 2 2 xik − yki + ηk2 ∇fπi (xik ) − ∇fπi (x⋆ ) + 2ηk Bπi (yki , x⋆ ) − Bπi (xik , x⋆ ) k k k k 2 xik − yki + 2 ηk2 Lπi − ηk Bπi (xik , x⋆ ) + 2ηk Bπi (yki , x⋆ )
xi+1 − yki+1 k = = (a)
≤ (b)
≤ (c)
2
k
≤ xik − yki
2
k
k
− ηk Bπi (xik , x⋆ ) + 2ηk Bπi (yki , x⋆ ), k
k
1 ⇒ 2L̂ k k η n+1 2 1 k ηk Lπi ≤ 2 in (c). Finally, we sum the above inequality from i = 1 to n, use xk+1 = xk and k ykn+1 = yk1 = x⋆ (see (9)), and rearrange terms to complete the proof.
where we use Bπi (yki , xik ) ≥ 0 in (a), apply Lemma 1 to fπi in (b), and notice that ηk ≤
E.2. Proof of Lemma 7 Proof Expanding the definitions of B and Bπi , we know k
B(xik , x⋆ ) − Bπi (xik , x⋆ ) =f (xik ) − f (x⋆ ) − ∇f (x⋆ ), xik − x⋆ k E D − fπi (xik ) − fπi (x⋆ ) − ∇fπi (x⋆ ), xik − x⋆ k
k
k
=f (xik ) − fπi (xik ) − f (x⋆ ) + fπi (x⋆ ) k k D E − ∇f (x⋆ ) − ∇fπi (x⋆ ), xik − x⋆ . k
h i h i Since E fπi (x⋆ ) = f (x⋆ ) and E ∇fπi (x⋆ ) = ∇f (x⋆ ) under RR, after taking expectations on k k both sides, we obtain h i h D Ei E B(xik , x⋆ ) − Bπi (xik , x⋆ ) = E f (xik ) − fπi (xik ) − ∇f (x⋆ ) − ∇fπi (x⋆ ), xik . (21) k
k
k
Now, we denote by ℓi (x) ≜ fi (x) − ⟨∇fi (x⋆ ), x⟩ , ∀i ∈ [n] and ℓ(x) ≜ n1 ⟨∇f (x⋆ ), x⟩. Then, (21) implies that,
Pn
i=1 ℓi (x) = f (x) −
i h i h i (a) 1 X h E B(xik , x⋆ ) − Bπi (xik , x⋆ ) = E ℓ(xik ) − ℓπi (xik ) = E ℓπi (xik (i, j)) − ℓπi (xik ) , k k k k n j<i
(22) where (a) holds by Lemma 2. 22
R ANDOM R ESHUFFLING D OMINATES S TOCHASTIC G RADIENT D ESCENT
Note that ℓπi is Lπi -smooth by its definition and Assumption 3, we therefore have, for any k k j ∈ [i − 1], D E Lπi 2 ℓπi (xik (i, j)) − ℓπi (xik ) ≤ ∇ℓπi (xik ), xik (i, j) − xik + k xik (i, j) − xik k k k 2 D E Lπi 2 = ∇fπi (xik ) − ∇fπi (x⋆ ), xik (i, j) − xik + k xik (i, j) − xik k k 2 ∇fπi (xik ) − ∇fπi (x⋆ )
(b)
k
≤
2 2
+ Lπi xik (i, j) − xik
k
2Lπi
k
k
(c)
≤ Bπi (xik , x⋆ ) + Lπi xik (i, j) − xik k
2
k
(d)
2
≤ Bπi (xik , x⋆ ) + ηk2 Lπi ∇fπi (xjk ) − ∇fπj (xjk ) k
k
k
,
(23)
k
where (b) is by Cauchy-Schwarz inequality and AM-GM inequality, (c) is due to Lemma 1, and (d) holds by Lemma 3. Furthermore, we can bound ∇fπi (xjk ) − ∇fπj (xjk ) k
≤4
2
k
2
∇fπi (xjk ) − ∇fπi (x⋆ ) k k
+ 4 ∇fπj (xjk ) − ∇fπj (x⋆ ) 2
k
k
k
k
2
+ 4 ∇fπi (x⋆ ) k
k
k
≤8Lπi Bπi (xjk , x⋆ ) + 8Lπj Bπj (xjk , x⋆ ) + 4 ∇fπi (x⋆ ) k
2
+ 4 ∇fπj (x⋆ )
+ 4 ∇fπj (x⋆ )
2
k
2
,
k
where the last step is by Lemma 1 again. Plug the above inequality back into (23) and use ηk ≤ 1 2L̂ to have ℓπi (xik (i, j)) − ℓπi (xik ) ≤Bπi (xik , x⋆ ) + 2Bπi (xjk , x⋆ ) + 2Bπj (xjk , x⋆ ) k
k
k
k
+ 4ηk2 Lπi ∇fπi (x⋆ ) k
k
k
2
+ 4ηk2 Lπi ∇fπj (x⋆ ) k
2
.
(24)
k
Combine (22) and (24) to obtain i 1X h h i E B(xik , x⋆ ) − Bπi (xik , x⋆ ) ≤ E Bπi (xik , x⋆ ) + 2Bπi (xjk , x⋆ ) + 2Bπj (xjk , x⋆ ) k k k k n j<i 2 2 4ηk2 X + . E Lπi ∇fπi (x⋆ ) + Lπi ∇fπj (x⋆ ) k k k k n j<i
h i h i When j < i, we observe that E Bπi (xjk , x⋆ ) = E Bπj (xjk , x⋆ ) and the following two inequalities k k hold Pn 2 Ll ∥∇fl (x⋆ )∥2 = l=1 ≤ L̂σ⋆2 ≤ nL̄σ⋆2 , E Lπi ∇fπi (x⋆ ) k k n Pn nL̄−Ll 2 2 l=1 n−1 ∥∇fl (x⋆ )∥ E Lπi ∇fπj (x⋆ ) ≤ 2L̄σ⋆2 . = k k n Put everything together, rearrange terms, and use 1 ≤ n to conclude the desired inequality.
23
L IU
E.3. Proof of Lemma 8 Proof By smoothness (i.e., Assumption 3), we have Bπi (yki , x⋆ ) ≤ k
Lπi
k
2
2 (8) Lπki yki − x⋆ =
2
2 2 (8) ηk Lπki yki − yk1 =
2
i−1 X j=1
2
∇fπj (x⋆ )
.
k
Therefore, we can bound
E
" n X i=1
2 n i−1 2 X X η η 2 n2 L̄σ⋆2 Bπi (yki , x⋆ ) ≤ k E , Lπi ∇fπj (x⋆ ) ≤ k k k k 2 6 #
i=1
j=1
where the last step is due to Lemma E.1 of Liu and Zhou (2024) (the constant here is slightly better, since ∇f (x⋆ ) = 0 in our setting leads to a provable improvement).
24