Stochastic Optimization
Part II: Faster online stochastic methods and Batch
stochastic optimization† ‡ Taiji Suzuki† Tokyo Institute of Technology
Graduate School of Information Science and Engineering
Department of Mathematical and Computing Sciences‡ JST, PRESTO

MLSS2015@Kyoto

1 / 59Outline

1

Getting stochastic gradient methods fasterBregman divergence and AdaGrad
Acceleration of stochastic gradient methods
Minimax optimality of first order online stochastic methods

2

Batch stochastic methods
Dual method: stochastic dual coordinate ascentPrimal method: SVRG, SAG and SAGA
Minimax optimality of first order batch stochastic methods

2 / 59Outline

1

Getting stochastic gradient methods fasterBregman divergence and AdaGrad
Acceleration of stochastic gradient methods
Minimax optimality of first order online stochastic methods

2

Batch stochastic methods
Dual method: stochastic dual coordinate ascentPrimal method: SVRG, SAG and SAGA
Minimax optimality of first order batch stochastic methods

3 / 59Outline

1

Getting stochastic gradient methods fasterBregman divergence and AdaGrad
Acceleration of stochastic gradient methods
Minimax optimality of first order online stochastic methods

2

Batch stochastic methods
Dual method: stochastic dual coordinate ascentPrimal method: SVRG, SAG and SAGA
Minimax optimality of first order batch stochastic methods

4 / 59Changing the metric (divergence)min L(x) + ψ(x)x{}1x (t) = argmin ⟨gt , x⟩ + ψ(x) + ∥x − x (t−1) ∥2
2η
x∈Rp

5 / 59Changing the metric (divergence)min L(x) + ψ(x)x{}1x (t) = argmin ⟨gt , x⟩ + ψ(x) + ∥x − x (t−1) ∥2Ht
2ηx∈Rp
∥x∥2H := x ⊤ Hx.

5 / 59Changing the metric (divergence)min L(x) + ψ(x)x{}1x (t) = argmin ⟨gt , x⟩ + ψ(x) + ∥x − x (t−1) ∥2Ht
2ηx∈Rp
∥x∥2H := x ⊤ Hx.
Choice of Ht
Hessian Ht = ∇∇⊤ L(x (t−1) ): Newton method
Fisher information matrix Ht = EZ |x (t−1) [−∇x ∇⊤x px (Z )|x=x (t−1) ]:Natural gradient(x is a parameter of a parametric model {px }x )c.f. Bregman divergence.
Bϕ (x||x ′ ) := ϕ(x) − ϕ(x ′ ) − ⟨∇ϕ(x ′ ), x − x ′ ⟩.
→ Mirror descent

5 / 596 / 59

AdaGrad (Duchi et al. (2011))
Let1

Ht = Gt2 + δI
for some δ ≥ 0, where Gt is either of the followings:(Full)Gt =t∑gτ gτ⊤ ,τ =1(Diag)

Gt = diag( t∑)gτ gτ⊤.

τ =1

AdaGrad stretches flat directions and shrinks steep directions.
Ada-SGD:{}
1x (t) = argmin ⟨gt , x⟩ + ψ(x) + ∥x − x (t−1) ∥2Ht .
2η
x∈Rp
∑
t
Ada-SRDA: for ḡt = 1t τ =1 gτ ,{}1∥x∥2Ht .
x (t) = argmin ⟨ḡt , x⟩ + ψ(x) +2tη
x

7 / 59Analysis of AdaGrad
Theorem
Let q = 2 for FULL, and q = ∞ for Diag. Define the regret as)
1 ∑(Q(T ) :=
ℓt+1 (x (t) ) + ψ(x (t) ) − ℓt (β ∗ ) − ψ(β ∗ ) .
T
T

t=1

Ada-SGD: ∀δ ≥ 0,
Q(T ) ≤

maxt≤T {∥x ∗ − x (t) ∥2q }/η + 2η [ 1/2 ]δ∥x ∗ ∥22 +tr GT .
Tη
2T

Ada-SRDA: for δ ≥ maxt ∥gt ∥2 ,Q(T ) ≤∥x ∗ ∥2q /η + 2η [ 1/2 ]δ∥β ∗ ∥22 +tr GT .
Tη2T
8 / 59Analysis of AdaGrad
SupposeThe gradient is unbalanced:|gt,j |2 ≤ Gj −2log(p)
E[L(x (T ) ) − L(x ∗ )] ≤ C √T(Ada-SGD)(ordinary SGD)

(j = 1, . . . , p, ∀t).

E[L(x

(T )

√
p
G E[maxt ∥x (t) ∥]√) − L(x )] ≤ C
≤C√
T
T

√
p

∗

→

log(p)

Much improvement.
9 / 59

AdaGrad is used in various applications including sparse learning and deeplearning.
In deep learning, we often encounter a phenomenon called plateau, that is,
we are stuck in a flat region.
It is hard to get out from plateau by standard SGD. AdaGrad adaptively
adjust the search space to get out of plateau.
AdaGrad is one of the standard optimization methods for deep learning.
Related methods: AdaDelta (Zeiler, 2012), RMSProp (Tieleman and
Hinton, 2012), Adam (Kingma and Ba, 2014).

10 / 59

Outline

1

Getting stochastic gradient methods fasterBregman divergence and AdaGrad
Acceleration of stochastic gradient methods
Minimax optimality of first order online stochastic methods

2

Batch stochastic methods
Dual method: stochastic dual coordinate ascentPrimal method: SVRG, SAG and SAGA
Minimax optimality of first order batch stochastic methods

11 / 59

Nesterov’s acceleration of SGD
Assumption:
the expected loss L(x) is γ-smooth.
the variance of gradient is bounded by σ 2 :EZ [∥∇β ℓ(Z , β) − ∇L(β)∥2 ] ≤ σ 2 .
→ combining with Nesterov’s acceleration, the convergence can be got
faster.
Acceleration for SGD: Hu et al. (2009)
Accleration for SRDA: Xiao (2010), Chen et al. (2012)
General method and analysis (including non-convex): Lan (2012),
Ghadimi and Lan (2012, 2013)
Ez1:T [Lψ (x (T ) )] − Lψ (x ∗ ) ≤ C

(

2
σD
√
+ DT 2γT)(D is the diameter: E[∥x (t) − x ∗ ∥2 ] ≤ D 2 (∀t))
12 / 59

Speed up of accelerated SGD
Ez1:T [Lψ (x

(T )

(

∗

)] − Lψ (x ) ≤ C

σD
D 2γ
√ + 2
T
T

)

σ 2 is the variance of the gradient estimate:EZ [∥∇β ℓ(Z , β) − ∇L(β)∥2 ] ≤ σ 2 .

The variance can be reduced by simply taking average:
1 ∑g=∇ℓ(zk , x (t−1) )KKg = ∇ℓ(z, x(Variance)

σ2

(t−1)

)

⇒

k=1

σ 2 /K

Computing independent gradients can be parallelized.
As σ → 0, the bound goes to O(1/T 2 ): non-stochastic Nesterov’s
acceleration.
13 / 59(a) Objective for L1(b) Objective for Elastic-net

Numerical comparison on synthetic data with (a) L1 regularization (Lasso)and (b) Elastic-net regularization (figure is from Chen et al. (2012)).
SAGE: Accelerated SGD (Hu et al., 2009), AC-RDA: Accelerated stochastic RDA (Xiao,2010), AC-SA: Accelerated stochastic approximation Ghadimi and Lan (2012), ORDA:
Optimal stochastic RDA (Chen et al., 2012)14 / 59

Accelerated SA for strongly convex objective
Assumption: Objective is µ-strongly convex and γ-smooth.
Accelerated stochastic approximation: Hu et al. (2009), Ghadimi andLan (2012)
Ez1:T [Lψ (x

(T )

(

∗

)] − Lψ (x ) ≤ C

σ2
γR 2
+ 2
µT
T

)
.

Multi-stage accelerated stochastic approximation: Chen et al. (2012),
Ghadimi and Lan (2013)
Ez1:T [Lψ (x

(T )

∗

)] − Lψ (x ) ≤ C

(

(
))
√
σ2
µ
+ exp −CTµT
γ

σ = 0 gives the batch optimal rate.
15 / 59

Summary of convergence rates
Online methods (expected risk minimization):GR√ (non-smooth, non-strongly convex) Polyak-Ruppert averaging
TG2(non-smooth, strongly convex)
Polynomial averaging
µT
σR
R 2L
√ + 2 (smooth, non-strongly convex)
Acceleration
T( √
T
)
σ2
µ
+ exp −T (smooth, strongly convex)
AccelerationµT
L

G : upper bound of norm of gradient, R: diameter of the domain,
L: smoothness, µ: strong convexity, σ: variance of the gradient

16 / 59

Outline

1

Getting stochastic gradient methods fasterBregman divergence and AdaGrad
Acceleration of stochastic gradient methods
Minimax optimality of first order online stochastic methods

2

Batch stochastic methods
Dual method: stochastic dual coordinate ascentPrimal method: SVRG, SAG and SAGA
Minimax optimality of first order batch stochastic methods

17 / 59

Minimax optimal rate of stochastic first order methods
min L(x) = min EZ [ℓ(Z , x)]
x∈B

x∈B

Conditionĝx ∈ ∂x ℓ(Z , x) is bounded as ∥E[ĝx ]∥ ≤ G (∀x ∈ B).
The domain B contains a ball with radius R.
L(x) is µ-strongly convex (µ = 0 is allowed).

Theorem (Minimax optimality (Agarwal et al., 2012, Nemirovsky and Yudin,1983))For any first order algorithm, there exist loss function ℓ and distribution
P(Z ) satisfying the assumption on which the algorithm must suﬀer{}
GR G 2 GR(T )∗E[L(x ) − L(x )] ≥ c min √ ,
,√
.
p
T µT
SGD and SRDA achieve this optimal rate.
First order algorithm: an algorithm that depends on only the loss and its
gradient (ℓ(Z , x), ĝx ) for a query point x. (SGD, SRDA are included.)
18 / 59

Outline

1

Getting stochastic gradient methods fasterBregman divergence and AdaGrad
Acceleration of stochastic gradient methods
Minimax optimality of first order online stochastic methods

2

Batch stochastic methods
Dual method: stochastic dual coordinate ascentPrimal method: SVRG, SAG and SAGA
Minimax optimality of first order batch stochastic methods

19 / 59

From expectation to finite sum∫

Online:P(x) = E[ℓ(Z , x)] =ℓ(Z , x)dP(Z )↓Batch:

1∑
P(x) =ℓ(zi , x)nni=120 / 59

From online to batch
In the batch setting, the data are fixed. We just minimize the objectivefunction defined by1∑ℓi (x) + ψ(x).
n i=1n

P(x) =


We construct a method thatuses few observations per iteration (like online method),converges linearly (unlike online method):
T > (n + γ/λ)log(1/ϵ)
to achieve ϵ accuracy for γ-smooth loss andλ-strongly convex regularization.


21 / 59

Three methods that must be remembered
Stochastic Average Gradient descent, SAG (Le Roux et al., 2012,
Schmidt et al., 2013, Defazio et al., 2014)

Stochastic Variance Reduced Gradient descent, SVRG (Johnson
and Zhang, 2013, Xiao and Zhang, 2014)

Stochastic Dual Coordinate Ascent, SDCA (Shalev-Shwartz andZhang, 2013a)- SAG and SVRG are methods performed on the primal.
- SDCA is on the dual.

22 / 59

Assumptions
P(x) = n1|∑ni=1 ℓi (x) +{z}smoothψ(x)| {z }
strongly convex

Assumption:ℓi : Loss is γ-smooth.
√ψ: reg func is λ-strongly convex. Typically λ = O(1/n) or O(1/ n).
Example:Loss function

Regularization functionsmoothed hinge loss

L2 regularization

logistic loss

Elastic net regularization
ψ̃(x) + λ∥x∥2 (with small h)0023 / 59

Outline

1

Getting stochastic gradient methods fasterBregman divergence and AdaGrad
Acceleration of stochastic gradient methods
Minimax optimality of first order online stochastic methods

2

Batch stochastic methods
Dual method: stochastic dual coordinate ascentPrimal method: SVRG, SAG and SAGA
Minimax optimality of first order batch stochastic methods

24 / 59

Coordinate Descent

25 / 59

Note on CD

failuresuccess

Left hand side: CD fails. No descent direction.
To make CD success, the∑objective should have descent direction.
Ideally, separable f (x) = pj=1 fj (xj ).

26 / 59

Coordinate descent in primalmin{P(x)} = min{f (x) + ψ(x)} = min{f (x) +
x

x

x

p
∑

ψj (xj )}j=1

Coordinate descent (sketch)
1 Choose j ∈ {1, . . . , p} in some way．(typically, random choice)
2 j-th coordinate xj is updated so that the objective is decreased,
Usually a block of coordinates are updated instead of one coordinate(block coordinate descent).

27 / 59

Coordinate descent in primalmin{P(x)} = min{f (x) + ψ(x)} = min{f (x) +
x

x

x

p
∑

ψj (xj )}j=1

Coordinate descent (sketch)
1 Choose j ∈ {1, . . . , p} in some way．(typically, random choice)
2 j-th coordinate xj is updated so that the objective is decreased, e.g.,
(t)

(t−1)

xj ← argminxj P(x1or(x (t) )for gj = ∂f ∂x
j
(t+1)

xj

(t−1)

, . . . , xj , . . . , xp),(t−1) 2← argminxj ⟨gj , xj ⟩ + ψj (xj ) + 2η1 t ∥xj − xj∥ .

Usually a block of coordinates are updated instead of one coordinate(block coordinate descent).
27 / 59

Convergence of primal CD method
We consider a separable regularization:minx {P(x)} = minx {f (x) + ψ(x)} = minx {f (x) +∑p

j=1 ψj (xj )}.

Assumption: f is γ-smooth (∥∇f (x) − ∇f (x ′ )∥ ≤ γ∥x − x ′ ∥)Cyclic (Saha and Tewari, 2013, Beck and Tetruashvili, 2013)
P(x (t) ) − R(x ∗ ) ≤γp∥x (0) − x ∗ ∥2
= O(1/t) (with isotonicity).
2t

Random choice (Nesterov, 2012, Richtárik and Takáč, 2014)
No acceleration: O(1/t).
Nesterov’s acceleration: O(1/t 2 ) (Fercoq and Richtárik, 2013).
f is α-strongly convex: O(exp(−C (α/γ)t)). √f is α-strongly conv + acceleration: O(exp(−C α/γt)) (Lin et al.,2014).

Nice review is given by Wright (2015).
28 / 59

Stochastic Dual Coordinate Ascent, SDCA
Suppose that ∃fi : R → R such that ℓ(zi , x) = fi (ai⊤ x).
Let A = [a1 , . . . , an ].
(Primal)


{
inf p

x∈R

1∑
fi (ai⊤ x) + ψ(x)
n
n

}



i=1



29 / 59

Stochastic Dual Coordinate Ascent, SDCA
Suppose that ∃fi : R → R such that ℓ(zi , x) = fi (ai⊤ x).
Let A = [a1 , . . . , an ].
{(Primal)

inf p

x∈R

1∑
fi (ai⊤ x) + ψ(x)
n
n

}



i=1





[Fenchel’s duality theorem]
inf {f (A⊤ x) + nψ(x)} = − inf n {f ∗ (y ) + nψ ∗ (−Ay /n)}
y ∈R

x∈Rp



{

(Dual)

inf n

y ∈R

1
n

n
∑
i=1

(
)}
1
fi ∗ (yi ) + ψ ∗ − Ay
n







We used the following facts:∑∑
For f (α) = ni=1 fi (αi ), we have f ∗ (β) = ni=1 fi ∗ (βi ).
For ψ̃(x) = nψ(x), we have ψ̃ ∗ (y ) = nψ ∗ (y /n).

29 / 59

Remarks
{
sup

y ∈Rn

)}
(
n
1∑ ∗
1
fi (yi ) + ψ ∗ − Ay
n
n
i=1

∑
The dual loss term ni=1 fi ∗ (yi ) is separable.
Each coordinate yi aﬀects the objective through only the i-th data:fi ∗ (yi ),)(1∗ψ − (a1 y1 + · · · + ai yi + · · · + an yn ) .
n→ Coordinate descent behaves like online methods!

The loss fi is smooth ⇔ fi ∗ is strongly convex.
The reg func ψ is strongly convex ⇔ ψ ∗ is smooth.
30 / 59

Algorithm of SDCA
SDCA (Shalev-Shwartz and Zhang, 2013a)
Iterate the following for t = 1, 2, . . .
1
2

Pick up an index i ∈ {1, . . . , n} uniformly at random.
Update the i-th coordinate yi so that the objective function isdecreased.

31 / 59

Algorithm of SDCA
SDCA (Shalev-Shwartz and Zhang, 2013a)
Iterate the following for t = 1, 2, . . .
1
2

Pick up an index i ∈ {1, . . . , n} uniformly at random.
Update the i-th coordinate yi :(let A\i = [a1 , . . . , ai−1 , ai+1 , . . . , an ], and y\i = (yj )j̸=i )
(t)

• yi

(t)

• yj

{
( 1
)
(t−1)
∈ argmin fi ∗ (yi ) + nψ ∗ − (ai yi + A\i y\i )
n
yi ∈R
}
1
(t−1) 2+ ∥yi − yi
∥ ,
2η
(t−1)

=yj

(for j ̸= i).

31 / 59

Algorithm of SDCA
SDCA (linearized version) (Shalev-Shwartz and Zhang, 2013a)
Iterate the following for t = 1, 2, . . .
1

Pick up an index i ∈ {1, . . . , n} uniformly at random.

2

Calculate x (t−1) = ∇ψ ∗ (−Ay (t−1) /n).

3

Update the i-th coordinate yi :
(t)

• yi

(t)

• yj

}
{
1
(t−1) 2∈ argmin fi ∗ (yi ) − ⟨x (t−1) , ai yi ⟩ + ∥yi − yi
∥
2η
yi ∈R
(t−1)

=yj

(for j ̸= i).

If the reg func ψ is λ-strongly covnex, ψ ∗ is 1/λ-smooth and thus
diﬀerentiable: x (t) = ∇ψ ∗ (−Ay (t) /n).
x (t) is actually the primal variable.
Computational complexity per iteration is same as online methods!
Important relation: prox(q|g ∗ ) = q − prox(q|g ). primal!

31 / 59

Algorithm of SDCA
SDCA (linearized version) (Shalev-Shwartz and Zhang, 2013a)
Iterate the following for t = 1, 2, . . .
1

Pick up an index i ∈ {1, . . . , n} uniformly at random.

2

Calculate x (t−1) = ∇ψ ∗ (−Ay (t−1) /n).

3

Update the i-th coordinate yi :
(t)

• yi

}
{
1
(t−1) 2∈ argmin fi ∗ (yi ) − ⟨x (t−1) , ai yi ⟩ + ∥yi − yi
∥
2η
yi ∈R
(t−1)

= prox(yi
(t)

• yj

(t−1)

=yj

+ ηai⊤ x (t−1) |ηfi ∗ ),(for j ̸= i).

If the reg func ψ is λ-strongly covnex, ψ ∗ is 1/λ-smooth and thus
diﬀerentiable: x (t) = ∇ψ ∗ (−Ay (t) /n).
x (t) is actually the primal variable.
Computational complexity per iteration is same as online methods!
Important relation: prox(q|g ∗ ) = q − prox(q|g ). primal!
31 / 59

Convergence analysis of SDCA
Assuption:fi is γ-smooth.
ψ is λ-strongly convex.

TheoremSuppose there exists R such that ∥ai ∥ ≤ R. Then, for η = λn/R 2 , we have
(
)
(
)
R 2γ
T
(T )
(T )
E[P(x̄ ) − D(ȳ )] ≤ n +exp −(D(y (0) ) − D(y ∗ )).
R 2γ
λ
n+
λ

E[·] is taken w.r.t. the choice of coordinates.

Linear convergence!
Required number of iterations to achieve ϵ:)(R 2γlog ((n + γ/λ)/ϵ) .
T ≥C n+λ32 / 59

Comparison with the non-stochastic method
How much computation is required to achieve E[P(x̄ (T ) ) − P(x ∗ )] ≤ ϵ?
Let κ = γ/λ (condition number).
SDCA:(n + κ) log ((n + κ)/ϵ)Ω((n + κ) log (1/ϵ)) iterations ×
Non-stochastic first order method:Ω(1) per iterationnκ log (1/ϵ)Ω(κ log (1/ϵ)) iterations

×Ω(n) per iteration

Sample size n = 100, 000, reg param λ = 1/1000, smoothness γ = 1:
n×κ = 108 ,n+κ = 105 .
33 / 59

Numerical comparison between SDCA, SDCA-perm (randomly shuﬄed
cyclic), SGD (figure is from Shalev-Shwartz and Zhang (2013a)).
34 / 59

Nesterov’s acceleration of SDCA
Accelerated SDCA (Lin et al., 2014)
Set α = n1√λγ.
(t−1) +αw (t−1)

1y (t) = ȳ

2

Pick up an index i ∈ {1, . . . , n} uniformly at random.

3

Calculate x (t−1) = ∇ψ ∗ (−Ay (t−1) /n).

41+αUpdate the i-th coordinate:
{
}
(t)
(t)
(t−1) 2• wi ∈ argminwi ∈R fi∗ (wi )−⟨x (t−1) ,ai wi ⟩+ αn∥w−y−(1−α)w
∥
i
i
i
2γ

5

(t)

(t−1)

(t)

=(1 − α)wj

(t)

(t)

• wj

(t)

(for j ̸= i).

+ yi

(t−1)

ȳi = yi + nα(wi − (1 − α)wi(t)(t)
ȳj = yj(for j ̸= i).

(t)− yi ),

Shalev-Shwartz and Zhang (2014) also proposed a double-loop acceleration

35 / 59

Convergence of accelerated SDCA
∥ai ∥ ≤ R (∀i) ∥A∥: spectral norm of A

TheoremConvergence of acc. SDCA If
(
√
T ≥

n+

γnR 2
λ

)

(
log

C γ∥A∥22
λnϵ

)
,

then
(Duality gap)

(normal)(accelerated)

E[P(x (T ) ) − D(y (T ) )] ≤ ϵ.
(γ)
n+log((n + κ)/ϵ)
λ
√ )
(
γn
n+
log((n + κ)/ϵ)λ
36 / 59

Mini-batch SDCA
Instead of choosing one coordinate yi , we may choose a block of
coordinates yI where I ⊆ {1, . . . , n}.
Typically, {1,∪. . . , n} is divided into K equally sized groups: I1 , . . . , IK s.t.
|Ik | = n/K , k Ik = {1, . . . , n}, Ik ∩ Ik ′ = ∅.

Mini-batch technique (Takáč et al., 2013, Shalev-Shwartz and Zhang, 2013b).
If K = n, we observe only one data at each iteration.
37 / 59

Mini-batch SDCA
Mini-batch SDCA (stochastic block coordinate descent)
For t = 1, 2, . . . ,, iterate the following:1

Randomly pick up a mini-batch I ⊆ {1, . . . , n} so thatP(i ∈ I ) = 1/K (∀i).

2x (t−1) = ∇ψ ∗ (−Ay (t−1) /n).

3

Update y (t) as
(t)

• yI

∈ argminyi (i∈I )

(t)
(t−1)
• yi =yi|I |{∑fi ∗ (yi ) − ⟨x (t−1) , AI yI ⟩ +

i=1

}
1
(t−1) 2∥yI − yI∥ ,2η(i ̸∈ I ).

The update of yi can be parallelized:(t−1)

yi = prox(yi+ ηai⊤ x (t−1) |ηfi ∗ ) (i ∈ I ).
38 / 59

Mini-batch SDCA
Mini-batch SDCA (stochastic block coordinate descent)
For t = 1, 2, . . . ,, iterate the following:1

Randomly pick up a mini-batch I ⊆ {1, . . . , n} so thatP(i ∈ I ) = 1/K (∀i).

2x (t−1) = ∇ψ ∗ (−Ay (t−1) /n).

3

Update y (t) as
(t)

• yI

∈ argminyi (i∈I )

(t)
(t−1)
• yi =yi

|I |
{∑
[

fi ∗ (yi ) − ⟨x (t−1) , Ai yi ⟩ +

i=1

}
1
(t−1) 2 ]∥yi − yi∥,2η(i ̸∈ I ).

The update of yi can be parallelized:(t−1)

yi = prox(yi+ ηai⊤ x (t−1) |ηfi ∗ ) (i ∈ I ).
38 / 59

Convergence of mini-batch SDCA
Assuption:fi is γ-smooth.
ψ is λ-strongly convex.

Theorem2
Suppose there exists R such that ∥A⊤
I AI ∥ ≤ R (∀I ). Then, for2
η = λn/R , we have
(
)
)
(
2γ
T
R
exp −
E[P(x̄ (T ) ) − D(ȳ (T ) )] ≤ K +(D(y (0) ) − D(y ∗ )).
2
λ
K+Rγ
λ

E[·] is taken w.r.t. the choice of coordinates.

(
)
R 2γ
T ≥C K+log((n + κ)/ϵ)λachieves ϵ accuracy. → iteration complexity is improved (if R 2 is not large
and parallelization is used).

39 / 59

Outline

1

Getting stochastic gradient methods fasterBregman divergence and AdaGrad
Acceleration of stochastic gradient methods
Minimax optimality of first order online stochastic methods

2

Batch stochastic methods
Dual method: stochastic dual coordinate ascentPrimal method: SVRG, SAG and SAGA
Minimax optimality of first order batch stochastic methods

40 / 59

Primal methods
The key idea: reduce the variance of gradient estimate.
1∑ℓi (x)ni=1| {z }n+ψ(x)How to approximate this?

41 / 59

Primal methods
The key idea: reduce the variance of gradient estimate.
1∑ℓi (x)ni=1| {z }n+ψ(x)How to approximate this?
⟨g ,x⟩

Online method: pick up î ∈ {1, . . . , n} randomly, and use linear
approximation.
∑g = ∇ℓî (x) ⇒ E[g ] = n1 ni=1 ∇ℓi (x)This is an unbiased estimator of the full gradient.
How about variance?
→ Variance is the problem!
→ In the batch setting, it is easy to reduce the variance.
41 / 59

Stochastic Variance Reduced Gradient descent,
SVRG (Johnson and Zhang, 2013, Xiao and Zhang, 2014)
minx {L(x) + ψ(x)} = minx { n1∑ni=1 ℓi (x) + ψ(x)}With fixed reference point x̂ which is close to x, a reduced variancegradient estimator is given as∑g = ∇ℓi (x)−∇ℓi (x̂) + n1 nj=1 ∇ℓj (x̂).
|{z}∇L(x̂)



Bias: unbiased,
E[g ] =

1∑1∑[∇ℓi (x) − ∇ℓi (x̂) + ∇L(x̂)] =∇ℓi (x) = ∇L(x).
n
n
n

n

i=1

i=1

Variance ?
42 / 59

A key observationg = ∇ℓi (x) − ∇ℓi (x̂) + ∇L(x̂).
Variance:∑n2i=1 ∥∇ℓi (x) − ∇ℓi (x̂) + ∇L(x̂) − ∇L(x)∥∑2= n1 ni=1 ∇ℓi (x) − ∇ℓi (x̂)∥2 − ∥∇L(x̂) − ∇L(x)22

Var[g ] = n1

1 ∑n

≤n

(∵ Var[X ] = E[∥X ∥ ] − ∥E[X ]∥ )2
i=1 ∥∇ℓi (x) − ∇ℓi (x̂)∥
2≤ γ∥x − x̂∥ .

The variance could be small if x and x̂ are close and ℓi is smooth.
Main strategy:
Calculate the full gradient at x̂.
Update xt several times, say, O(n) times.
Set x̂ = xt .
43 / 59

}

Variance

44 / 59

}

Variance}}

44 / 59

Algorithm of SVRG
The algorithm consists of inner loop and outer loop.

SVRG
For t = 1, 2, . . . , iterate the following:1

Set x̂ = x̂ (t−1) , x[0] = x̂,ĝ = ∇L(x̂) = n12i=1 ∇ℓi (x̂).

(full gradient)

For k = 1, . . . , m, execute the following:
1
2

3

3

∑n

Uniformly sample i ∈ {1, . . . , n}．
Setg = ∇ℓi (x[k−1] ) − ∇ℓi (x̂) + ĝ .
Update x[k] as
1 ∑mSet x̂ (t) = m(variance reduction)x[k] = prox(x[k−1] − ηg |ηψ).

k=1 x[k] .

Computational complexity until t iteration: O(t × (n + m)).
45 / 59

Convergence analysis
Assuption: ℓi is γ-smooth, and ψ is λ-strongly convex.

TheoremIf η and m satisfy η > 4γ and
4γ(m+1)ηρ := λ(1−4γ/η)m+ η(1−4γ/η)m< 1,

then, after T iteration the objective is bounded by
E[P(x̂ (T ) ) − P(x ∗ )] ≤ ρT (P(x̂ (0) ) − P(x ∗ ))The assumption of the theorem is satisfied by( )m ≥ Ω λγ .
Inner loop computation O(n + m) for each t.
Outer loop iteration T = O(log(1/ϵ)) until ϵ accuracy.
⇒ The whole computation :()γ
O ((n + m) log(1/ϵ)) = O (n + ) log(1/ϵ)λ

46 / 59

(c) rcv1(d) covertype(e) sido0

Numerical comparison between several stochastic methods on a batch
setting (figure is from Xiao and Zhang (2014)).
47 / 59

Related method: SAGASAGA (Defazio et al., 2014) does not require the double-loop,
but requires more memory.
Diﬀerence: the gradient estimate g
1∑∇ℓj (x̂j )nn(SAGA)

g = ∇ℓi (x (t−1) ) − ∇ℓi (x̂i ) +j=1(SVRG)

g = ∇ℓi (x (t−1) ) − ∇ℓi (x̂) +

n
1∑

n

∇ℓi (x̂)

j=1

x̂ depends on the data index i ∈ {1, . . . , n}.
x̂i is updated at every iteration:{(if i is chosen at the t-th round),x̂i = x (t−1)
x̂j is not changed (∀j ̸= i).
- Update rule of x (t) is same: x (t) = prox(x (t−1) − ηg |ηψ).
- We need to store all gradients ∇ℓi (x̂i ) (i = 1, . . . , n).
48 / 59

Algorithm of SAGASAGA (Defazio et al., 2014)
1

Pick up i ∈ {1, . . . , n} uniformly at random.

2

Update gj(t)(j = 1, . . . , n) as
{
(t)
gj =

3

∇ℓi (x (t−1) )(i = j),(t−1)
gj(otherwise).

Update x (t) as
1 ∑ (t−1)
gj
,
n
n

(t)

vt = gi

(t−1)

− gi

+

j=1

x

(t)

= prox(x(t−1)− ηvt |ηψ).

49 / 59

Convergence of SAGA
Assumption: ℓi is γ-smooth and λ-strongly convex (λ = 0 is allowed)(∀i = 1, . . . , n)．

TheoremSet η = 1/3γ. Thenλ = 0: for x̄ (T ) = T1

∑T

t=1 x

(t) , it holds that

E[P(x̄ (T ) ) − P(x ∗ )] ≤

4n
C0 .
T

λ > 0:
E[∥x

(T )

(
{
})T
1 λ
− x ∥ ] ≤ 1 − min
,
C1 .
4n 3γ
∗ 2

SAGA is adaptive to the strong convexity λ.
50 / 59

Stochastic Average Gradient descent, SAG
Like SAGA, SAG is also a single-loop method (Le Roux et al., 2012,
Schmidt et al., 2013).
Historically SAG was proposed earlier than SVRG and SAGA.
Proximal technique can not be involved in SAG→ SAGA was proposed to overcome this drawback.
1∑
P(x) =ℓi (x) ≃ ⟨g , x⟩.
n
n

i=1

∇ℓi (x (t−1) ) − ∇ℓi (x̂i ) 1 ∑+∇ℓj (x̂j )
g=
n
n
n

(SAG)

j=1

g is biased.
(SAGA)

g = ∇ℓi (x (t−1) ) − ∇ℓi (x̂i ) + n1

(SVRG)

g = ∇ℓi (x (t−1) ) − ∇ℓi (x̂) + n1∑n

j=1 ∇ℓj (x̂j )∑n

j=1 ∇ℓi (x̂)
51 / 59

Algorithm of SAG
SAG(0)

Initialize gi = 0 (i = 1, . . . , n).
For t = 1, 2, . . . , iterate the following:
1
2

3

Pick up i ∈ {1, . . . , n} uniformly at random.
Update gi ′ (i ′ = 1, . . . , n) as{∇ℓi (x (t−1) )(t)gi ′ =
(t−1)
gi ′
(t)

(i = i ′ ),(otherwise).

Update x (t) asη ∑ (t)gj .
nnx (t) = x (t−1) −j=1

52 / 59

Convergence analysis of SAG∑
Assumption: ℓi is γ-smooth and P(x) = n1 ni=1 ℓi (x) is λ-strongly
convex (λ = 0 is allowed)．
Milder condition than SAGA because the strong convexity is about P(x)rather than the loss function ℓi .

Theorem (Convergence rate of SAG)
1
Set η = 16γ. Then SAG converges as∑(t) に対し，
λ = 0: x̄ (T ) = T1 T
t=1 x

E[P(x̄ (T ) ) − P(x ∗ )] ≤

32n
C0
T

λ > 0:
E[∥x

(T )

})T
(
{
1 λ
C0 .
− x ∥ ] ≤ 1 − min,
8n 16γSAG also has adaptivity.

∗ 2

53 / 59

Catalyst: Acceleration of SVRG, SAG, SAGA
Catalyst (Lin et al., 2015)
Iterate the following for t = 1, 2, . . . :1

2

Find an approximated solution of a modified problem which has
higher strong convexity:{}
x (t) ≃ argminx P(x) + α2 ∥x − y (t−1) ∥2(up to ϵt precision).
Accelerate the solution: y (t) = x (t) + βt (x (t) − x (t−1) ).

− Catalyst is an acceleration
method of an inexact proximal point alg.
√
− For ϵt = C (1 − λ/2(λ + α))t ,
(
)t
√
λ
(t)
∗
′
P(x ) − P(x ) ≤ C 1 − 2(λ+α) .
− Using SVRG, SAG, SAGA with α = max{c γn , λ} in the inner loop√
achieves (n + γnλ ) log(1/ϵ) overall computation.
− This is a universal method but is sensitive to the choice of the inner
54 / 59
loop iteration number and α.

Summary and comparison of batch methods
Properties of the batch methods
Method
SDCA
SVRG
P/D
Dual
PrimalMemory eﬃciency✓✓
Acceleration (µ > 0)✓
Catalyst
Other remarkℓi (β) = fi (xi⊤ β) double loop

SAG
Primal△Catalyst
smooth reg.

55 / 59

Summary and comparison of batch methods
Properties of the batch methods
Method
SDCA
SVRG
P/D
Dual
PrimalMemory eﬃciency✓✓
Acceleration (µ > 0)✓
Catalyst
Other remarkℓi (β) = fi (xi⊤ β) double loop

MethodSAG
Primal△Catalyst
smooth reg.

Convergence rate (up to log term of γ, µ)λ>0
λ = 0 Acceleration(µ > 0)√SDCA(n + λγ ) log(1/ϵ)

-

(n +

) log(1/ϵ)

SVRG(n + λγ ) log(1/ϵ)

√λ

nγ

-

(n +

) log(1/ϵ)

SAG(n + λγ ) log(1/ϵ)

√λ

γn/ϵ

(n +

nγ

nγ
λ ) log(1/ϵ): Catalyst.
√
As for µ = 0, Catalyst gives an acceleration with convergence rate O(n γϵ ).
55 / 59

Outline

1

Getting stochastic gradient methods fasterBregman divergence and AdaGrad
Acceleration of stochastic gradient methods
Minimax optimality of first order online stochastic methods

2

Batch stochastic methods
Dual method: stochastic dual coordinate ascentPrimal method: SVRG, SAG and SAGA
Minimax optimality of first order batch stochastic methods

56 / 59

Minimax optimal convergence rate
Let κ = λγ be the condition number.
The iteration number
T ≥ (n + κ) log(1/ϵ)
is almost minimax, but not minimax.
The accelerated version(√ )T ≥ n + nκ log(1/ϵ)
is minimax up to log(1/ϵ) (Agarwal and Bottou, 2015).

57 / 59

Minimax optimality in the batch setting1∑
λℓi (x) + ∥x∥2n2nP(x) =

i=1

Assumption: ℓi is (γ − λ)-smooth (γ > λ).
First order oracle: for an input (x, i), it returns the pair (ℓi (x), ∇ℓi (x)).
First order algorithm: an algorithm that depends on only the return of thefirst order oracle for a query point x.
(SAG, SAGA, SVRG are included. SDCA is not included.)

Theorem (Minimax optimal rate for FOA (Agarwal and Bottou, 2015))
For any first order algorithm, there exist functions ℓi (i = 1, . . . , n)satisfying the assumption on which the algorithm must perform at least√T ≥ Ω(n + n(κ − 1) log(1/ϵ))
calls for the first order oracle to get ∥x (T ) − x ∗ ∥ ≤ ϵ∥x ∗ ∥.
58 / 59

Summary of part II
Getting faster the stochastic optimization methods
AdaGrad: adaptive choice of the step size
Acceleration of stochastic approximation2

σR
Non strongly conv: √+ RT 2LT( √)2σStrongly conv: µT+ exp − µL T

Batch optimization methods
SDCA, SVRG, SAG, SAGA
Variance
reduction()
exp − n+1 µ T convergenceL()
1
exp − √T
convergence with acceleration
nµ
n+

L

Next part: stochastic alternating direction method of multipliers forstructured regularization and distributed optimization
59 / 59

A. Agarwal and L. Bottou. A lower bound for the optimization of finite
sums. In the 32nd International Conference on Machine Learning, pages
78–86, 2015.
A. Agarwal, P. L. Bartlett, P. Ravikumar, and M. J. Wainwright.
Information-theoretic lower bounds on the oracle complexity ofstochastic convex optimization. IEEE Transcations on Information
Theory, 58(5):3235–3249, 2012.
A. Beck and L. Tetruashvili. On the convergence of block coordinate
descent type methods. SIAM Journal on Optimization, 23(4):2037–2060, 2013.
X. Chen, Q. Lin, and J. Pena. Optimal regularized dual averaging methods
for stochastic optimization. In F. Pereira, C. Burges, L. Bottou, and
K. Weinberger, editors, Advances in Neural Information Processing
Systems 25, pages 395–403. Curran Associates, Inc., 2012. URL
http://papers.nips.cc/paper/4543-optimal-regularized-dual-averaging-methods-for-stochas
pdf.
A. Defazio, F. Bach, and S. Lacoste-Julien. SAGA: A fast incremental
gradient method with support for non-strongly convex composite
59 / 59

objectives. In Z. Ghahramani, M. Welling, C. Cortes, N. Lawrence, and
K. Weinberger, editors, Advances in Neural Information Processing
Systems 27, pages 1646–1654. Curran Associates, Inc., 2014. URL
http://papers.nips.cc/paper/
5258-saga-a-fast-incremental-gradient-method-with-support-f
pdf.
J. Duchi, E. Hazan, and Y. Singer. Adaptive subgradient methods foronline learning and stochastic optimization. Journal of Machine
Learning Research, 12:2121–2159, 2011.
O. Fercoq and P. Richtárik. Accelerated, parallel and proximal coordinate
descent. Technical report, 2013. arXiv:1312.5799.
S. Ghadimi and G. Lan. Optimal stochastic approximation algorithms forstrongly convex stochastic composite optimization I: A genericalgorithmic framework. SIAM Journal on Optimization, 22(4):1469–1492, 2012.
S. Ghadimi and G. Lan. Optimal stochastic approximation algorithms forstrongly convex stochastic composite optimization, II: shrinking
procedures and optimal algorithms. SIAM Journal on Optimization, 23(4):2061–2089, 2013.
59 / 59

C. Hu, W. Pan, and J. T. Kwok. Accelerated gradient methods forstochastic optimization and online learning. In Y. Bengio,
D. Schuurmans, J. Laﬀerty, C. Williams, and A. Culotta, editors,
Advances in Neural Information Processing Systems 22, pages 781–789.
Curran Associates, Inc., 2009. URL http://papers.nips.cc/paper/
3817-accelerated-gradient-methods-for-stochastic-optimizatio
pdf.

R. Johnson and T. Zhang. Accelerating stochastic gradient descent using
predictive variance reduction. In C. Burges, L. Bottou, M. Welling,
Z. Ghahramani, and K. Weinberger, editors, Advances in Neural
Information Processing Systems 26, pages 315–323. Curran Associates,
Inc., 2013. URL http://papers.nips.cc/paper/
4937-accelerating-stochastic-gradient-descent-using-predicti
pdf.
D. Kingma and J. Ba. Adam: A method for stochastic optimization. arXiv
preprint arXiv:1412.6980, 2014.
G. Lan. An optimal method for stochastic composite optimization.
Mathematical Programming, 133(1-2):365–397, 2012.
59 / 59

N. Le Roux, M. Schmidt, and F. R. Bach. A stochastic gradient methodwith an exponential convergence rate for finite training sets. In
F. Pereira, C. Burges, L. Bottou, and K. Weinberger, editors, Advancesin Neural Information Processing Systems 25, pages 2663–2671. Curran
Associates, Inc., 2012. URL http://papers.nips.cc/paper/
4633-a-stochastic-gradient-method-with-an-exponential-conve
rate-for-finite-training-sets.pdf.
H. Lin, J. Mairal, and Z. Harchaoui. A universal catalyst for first-order
optimization. Technical report, 2015. arXiv:1506.02186.
Q. Lin, Z. Lu, and L. Xiao. An accelerated proximal coordinate gradientmethod and its application to regularized empirical risk minimization.
Technical report, 2014. arXiv:1407.1296.
A. Nemirovsky and D. Yudin. Problem complexity and method eﬃciency
in optimization. John Wiley, New York, 1983.
Y. Nesterov. Eﬃciency of coordinate descent methods on huge-scale
optimization problems. SIAM Journal on Optimization, 22(2):341–362,
2012.
P. Richtárik and M. Takáč. Iteration complexity of randomized
59 / 59

block-coordinate descent methods for minimizing a composite function.
Mathematical Programming, 144:1–38, 2014.
A. Saha and A. Tewari. On the non-asymptotic convergence of cyclic
coordinate descent methods. SIAM Journal on Optimization, 23(1):576–601, 2013.
M. Schmidt, N. Le Roux, and F. R. Bach. Minimizing finite sums with thestochastic average gradient, 2013. hal-00860051.
S. Shalev-Shwartz and T. Zhang. Stochastic dual coordinate ascentmethods for regularized loss minimization. Journal of Machine Learning
Research, 14:567–599, 2013a.
S. Shalev-Shwartz and T. Zhang. Accelerated mini-batch stochastic dual
coordinate ascent. In Advances in Neural Information Processing
Systems 26, 2013b.
S. Shalev-Shwartz and T. Zhang. Accelerated proximal stochastic dual
coordinate ascent for regularized loss minimization. In Proceedings of
The 31st International Conference on Machine Learning, pages 64–72,
2014.
M. Takáč, A. Bijral, P. Richtárik, and N. Srebro. Mini-batch primal and
59 / 59

dual methods for SVMs. In Proceedings of the 30th International
Conference on Machine Learning, 2013.
T. Tieleman and G. Hinton. Lecture 6.5-rmsprop: Divide the gradient by arunning average of its recent magnitude. COURSERA: Neural Networks
for Machine Learning, 2012.
S. J. Wright. Coordinate descent algorithms. Mathematical Programming,
151(1):3–34, 2015.
L. Xiao. Dual averaging methods for regularized stochastic learning andonline optimization. Journal of Machine Learning Research, 11:2543–2596, 2010.
L. Xiao and T. Zhang. A proximal stochastic gradient method with
progressive variance reduction. SIAM Journal on Optimization, 24:
2057–2075, 2014.
M. D. Zeiler. ADADELTA: an adaptive learning rate method. CoRR,
abs/1212.5701, 2012. URL http://arxiv.org/abs/1212.5701.

59 / 59