Training Neural Networks with Stochastic Hessian-Free Optimization
Authors: Ryan Kiros
ICLR 2013

Abstract: Hessian-free (HF) optimization has been successfully used for training deep autoencoders and recurrent networks. HF uses the conjugate gradient algorithm to construct update directions through curvature-vector products that can be computed on the same order of time as gradients. In this paper we exploit this property and study stochastic HF with small gradient and curvature mini-batches independent of the dataset size for classification. We modify Martens' HF for this setting and integrate dropout, a method for preventing co-adaptation of feature detectors, to guard against overfitting. On classification tasks, stochastic HF achieves accelerated training and competitive results in comparison with dropout SGD without the need to tune learning rates.

arXiv:1301.3641v3  [cs.LG]  1 May 2013
T raining Neural Networks with Stochastic
Hessian-Free Optimization
Ryan Kiros
Department of Computing Science
University of Alberta
Edmonton, AB, Canada
rkiros@ualberta.ca
Abstract
Hessian-free (HF) optimization has been successfully used for training deep au-
toencoders and recurrent networks. HF uses the conjugate gr adient algorithm to
construct update directions through curvature-vector pro ducts that can be com-
puted on the same order of time as gradients. In this paper we exploit this property
and study stochastic HF with gradient and curvature mini-batches independent of
the dataset size. We modify Martens’ HF for these settings and integrate dropout,
a method for preventing co-adaptation of feature detectors, to guard against over-
ﬁtting. Stochastic Hessian-free optimization gives an int ermediary between SGD
and HF that achieves competitive performance on both classi ﬁcation and deep
autoencoder experiments.
1 Introduction
Stochastic gradient descent (SGD) has become the most popul ar algorithm for training neural net-
works. Not only is SGD simple to implement but its noisy updat es often leads to solutions that are
well-adapt to generalization on held-out data [1]. Furthermore, SGD operates on small mini-batches
potentially allowing for scalable training on large datase ts. For training deep networks, SGD can
be used for ﬁne-tuning after layerwise pre-training [2] whi ch overcomes many of the difﬁculties
of training deep networks. Additionally, SGD can be augment ed with dropout [3] as a means of
preventing overﬁtting.
There has been recent interest in second-order methods for t raining deep networks, partially due
to the successful adaptation of Hessian-free (HF) by [4], an instance of the more general family
of truncated Newton methods. Second-order methods operate in batch settings with less but more
substantial weight updates. Furthermore, computing gradi ents and curvature information on large
batches can easily be distributed across several machines.Martens’ HF was able to successfully train
deep autoencoders without the use of pre-training and was later used for solving several pathological
tasks in recurrent networks [5].
HF iteratively proposes update directions using the conjug ate gradient algorithm, requiring only
curvature-vector products and not an explicit computationof the curvature matrix. Curvature-vector
products can be computed on the same order of time as it takes t o compute gradients with an addi-
tional forward and backward pass through the function’s com putational graph [6, 7]. In this paper
we exploit this property and introduce stochastic Hessian-free optimization (SHF), a variation of HF
that operates on gradient and curvature mini-batches indep endent of the dataset size. Our goal in
developing SHF is to combine the generalization advantagesof SGD with second-order information
from HF. SHF can adapt its behaviour through the choice of bat ch size and number of conjugate
gradient iterations, for which its behaviour either become s more characteristic of SGD or HF. Ad-
ditionally we integrate dropout, as a means of preventing co -adaptation of feature detectors. We
1
perform experimental evaluation on both classiﬁcation and deep autoencoder tasks. For classiﬁca-
tion, dropout SHF is competitive with dropout SGD on all task s considered while for autoencoders
SHF performs comparably to HF and momentum-based methods. M oreover, no tuning of learning
rates needs to be done.
2 Related work
Much research has been investigated into developing adaptive learning rates or incorporating second-
order information into SGD. [8] proposed augmenting SGD wit h a diagonal approximation of the
Hessian while Adagrad [9] uses a global learning rate while d ividing by the norm of previous gra-
dients in its update. SGD with Adagrad was shown to be beneﬁci al in training deep distributed
networks for speech and object recognition [10]. To complet ely avoid tuning learning rates, [11]
considered computing rates as to minimize estimates of the e xpectation of the loss at any one time.
[12] proposed SGD-QN for incorporating a quasi-Newton appr oximation to the Hessian into SGD
and used this to win one of the 2008 PASCAL large scale learning challenge tracks. Recently, [13]
provided a relationship between HF, Krylov subspace descent and natural gradient due to their use
of the Gauss-Newton curvature matrix. Furthermore, [13] ar gue that natural gradient is robust to
overﬁtting as well as the order of the training samples. Othe r methods incorporating the natural
gradient such as TONGA [14] have also showed promise on speeding up neural network training.
Analyzing the difﬁculty of training deep networks was done b y [15], proposing a weight initial-
ization that demonstrates faster convergence. More recent ly, [16] argue that large neural networks
waste capacity in the sense that adding additional units fail to reduce underﬁtting on large datasets.
The authors hypothesize the SGD is the culprit and suggest exploration with stochastic natural gra-
dient or stochastic second-order methods. Such results fur ther motivate our development of SHF.
[17] show that with careful attention to the parameter initi alization and momentum schedule, ﬁrst-
order methods can be competitive with HF for training deep au toencoders and recurrent networks.
We compare against these methods in our autoencoder evaluation.
Related to our work is that of [18], who proposes a dynamic adj ustment of gradient and curvature
mini-batches for HF with convex losses based on variance est imations. Unlike our work, the batch
sizes used are dynamic with a ﬁxed ratio and are initialized a s a function of the dataset size. Other
work on using second-order methods for neural networks include [19] who proposed using the Jacobi
pre-conditioner for HF, [20] using HF to generate text in recurrent networks and [21] who explored
training with Krylov subspace descent (KSD). Unlike HF, KSD could be used with Hessian-vector
products but requires additional memory to store a basis for the Krylov subspace. L-BFGS has
also been successfully used in ﬁne-tuning pre-trained deep autoencoders, convolutional networks
[22] and training deep distributed networks [10]. Other dev elopments and detailed discussion of
gradient-based methods for neural networks is described in [23].
3 Hessian-free optimization
In this section we review Hessian-free optimization, large ly following the implementation of
Martens [4]. We refer the reader to [24] for detailed development and tips for using HF.
We consider unconstrained minimization of a function f : Rn →R with respect to parameters θ.
More speciﬁcally, we assume f can be written as a composition f(θ) = L(F (θ)) where L is a
convex loss function and F (θ) is the output of a neural network with ℓ non-input layers. We will
mostly focus on the case when f is non-convex. Typically L is chosen to be a matching loss to
a corresponding transfer function p(z) = p(F (θ)). For a single input, the (i + 1) -th layer of the
network is expressed as
yi+1 = si(Wi yi + bi) (1)
where si is a transfer function, Wi is the weights connecting layersi and i + 1 and bi is a bias vector.
Common transfer functions include the sigmoid si(x) = (1 + exp(−x))−1, the hyperbolic tangent
si(x) = tanh(x) and rectiﬁed linear units si(x) = max(x, 0). In the case of classiﬁcation tasks, the
2
loss function used is the generalized cross entropy and softmax transfer
L(p(z), t) = −
k∑
j=1
tj log(p(zj)), p (zj) = exp(zj )/
k∑
l=1
exp(zl) (2)
where k is the number of classes, t is a target vector and zj the j-th component of output vector z.
Consider a local quadratic approximation Mθ(δ) of f around θ:
f(θ + δ) ≈Mθ(δ) = f(θ) + ∇f(θ)T δ + 1
2δT Bδ (3)
where ∇f(θ) is the gradient of f and B is the Hessian or an approximation to the Hessian. If f
was convex, then B ⪰ 0 and equation 3 exhibits a minimum δ∗. In Newton’s method, θk+1, the
parameters at iteration k + 1, are updated as θk+1 = θk + αkδ∗
k where αk ∈[0, 1] is the rate and δ∗
k
is computed as
δ∗
k = −B−1∇f(θk−1) (4)
for which calculation requires O(n3) time and thus often prohibitive. Hessian-free optimizatio n
alleviates this by using the conjugate gradient (CG) algorithm to compute an approximate minimizer
δk. Speciﬁcally, CG minimizes the quadratic objective q(δ) given by
q(δ) = 1
2 δT Bδ + ∇f(θk−1)T δ (5)
for which the corresponding minimizer of q(δ) is −B−1∇f(θk−1). The motivation for using CG
is as follows: while computing B is expensive, compute the product Bv for some vector v can be
computed on the same order of time as it takes to compute ∇f(θk−1) using the R-operator [6].
Thus CG can efﬁciently compute an iterative solution to the l inear system Bδk = −∇(f(θk−1))
corresponding to a new update direction δk.
When f is non-convex, the Hessian may not be positive semi-deﬁniteand thus equation 3 no longer
has a well deﬁned minimum. Following Martens, we instead use the generalized Gauss-newton
matrix deﬁned as B = JT L
′′
J where J is the Jacobian of f and L
′′
is the Hessian of L 1. So long
as f(θ) = L(F (θ)) for convex L then B ⪰ 0. Given a vector v, the product Bv = JT L
′′
Jv is
computed successively by ﬁrst computing Jv , then L
′′
(Jv ) and ﬁnally JT (L
′′
Jv ) [7]. To compute
Jv , we utilize the R-operator. The R-operator of F (θ) with respect to v is deﬁned as
Rv{F (θ)}= lim
ǫ→0
F (θ + ǫv) −F (θ)
ǫ = Jv (6)
Computing Rv{F (θ)}in a neural network is easily done using a forward pass by computing Rv{yi}
for each layer output yi. More speciﬁcally,
Rv{yi+1}= Rv{Wi yi + bi}s′
i = ( v(Wi)yi + v(bi) + WiR{yi})s′
i (7)
where v(Wi) is the components of v corresponding to parameters between layers i and i + 1 and
R{y1}= 0 (where y1 is the input data). In order to compute JT (L
′′
Jv ), we simply apply back-
propagation but using the vector L
′′
Jv instead of ∇L as is usually done to compute ∇f. Thus, Bv
may be computed through a forward and backward pass in the same sense that L and ∇f = JT ∇L
are.
As opposed to minimizing equation3, Martens instead uses an additional damping parameterλ with
damped quadratic approximation
ˆMθ(δ) = f(θ) + ∇f(θ)T δ + 1
2δT ˆBδ = f(θ) + ∇f(θ)T δ + 1
2δT (B + λI)δ (8)
Damping the quadratic through λ gives a measure of how conservative the quadratic approximation
is. A large value of λ is more conservative and as λ →∞ updates become similar to stochastic
gradient descent. Alternatively, a small λ allows for more substantial parameter updates especially
1While an abuse of deﬁnition, we still refer to “curvature-ve ctor products” and “curvature batches” even
when B is used.
3
along low curvature directions. Martens dynamically adjusts λ at each iteration using a Levenberg-
Marquardt style update based on computing the reduction ratio
ρ = ( f(θ + δ) −f(θ))/(Mθ(δ) −Mθ(0)) (9)
If ρ is sufﬁciently small or negative, λ is increased while if ρ is large then λ is decreased. The
number of CG iterations used to compute δ has a dramatic effect on ρ which is further discussed in
section 4.1.
To accelerate CG, Martens makes use of the diagonal pre-conditioner
P =
[
diag
( m∑
j=1
∇f(j)(θ) ⊙∇f(j)(θ)
)
+ λI
]ξ
(10)
where f(j)(θ) is the value off for datapoint j and ⊙denotes component-wise multiplication. P can
be easily computed on the same backward pass as computing ∇f.
Finally, two backtracking methods are used: one after optim izing CG to select δ and the other a
backtracking linesearch to compute the rate α. Both these methods operate in the standard way,
backtracking through proposals until the objective no longer decreases.
4 Stochastic Hessian-free optimization
Martens’ implementation utilizes the full dataset for computing objective values and gradients, and
mini-batches for computing curvature-vector products. Na ively setting both batch sizes to be small
causes several problems. In this section we describe these p roblems and our contributions in modi-
fying Martens’ original algorithm to this setting.
4.1 Short CG runs, δ-momentum and use of mini-batches
The CG termination criteria used by Martens is based on a measure of relative progress in optimizing
ˆMθ. Speciﬁcally, if xj is the solution at CG iteration j, then training is terminated when
ˆMθ(xj ) − ˆMθ(xj−k )
ˆMθ(xj )
< ǫ (11)
where k =max(10, j/10) and ǫ is a small positive constant. The effect of this stopping cri teria
has a dependency on the strength of the damping parameter λ, among other attributes such as the
current parameter settings. For sufﬁciently large λ, CG only requires 10-20 iterations when a pre-
conditioner is used. As λ decreases, more iterations are required to account for path ological curva-
ture that can occur in optimizing f and thus leads to more expensive CG iterations. Such behavio r
would be undesirable in a stochastic setting where preferen ce would be put towards having equal
length CG iterations throughout training. To account for this, we ﬁx the number of CG iterations to
be only 3-5 across training for classiﬁcation and 25-50 for training deep autoencoders. Let ζ denote
this cut-off. Setting a limit on the number of CG iterations i s used by [4] and [20] and also has a
damping effect, since the objective function and quadraticapproximation will tend to diverge as CG
iterations increase [24]. We note that due to the shorter num ber of CG runs, the iterates from each
solution are used during the CG backtracking step.
A contributor to the success of Martens’ HF is the use of infor mation sharing across iterations.
At iteration k, CG is initialized to be the previous solution of CG from iter ation k −1, with a
small decay. For the rest of this work, we denote this as δ-momentum. δ-momentum helps correct
proposed update directions when the quadratic approximati on varies across iterations, in the same
sense that momentum is used to share gradients. This momentu m interpretation was ﬁrst suggested
by [24] in the context of adapting HF to a setting with short CG runs. Unfortunately, the use of δ-
momentum becomes challenging when short CG runs are used. Given a non-zero CG initialization,
ˆMθ may be more likely to remain positive after terminating CG and assuming f(θ + δ) −f(θ) < 0,
means that the reduction ratio will be negative and thus λ will be increased to compensate. While
this is not necessarily unwanted behavior, having this occu r too frequently will push SHF to be
too conservative and possibly result in the backtracking li nesearch to reject proposed updates. Our
4
solution is to utilize a schedule on the amount of decay used o n the CG starting solution. This is
motivated by [24] suggesting more attention on the CG decay in the setting of using short CG runs.
Speciﬁcally, if δ0
k is the initial solution to CG at iteration k, then
δ0
k = γeδζ
k−1, γ e = min(1.01γe−1, .99) (12)
where γe is the decay at epoch e, δ0
1 = 0 and γ1 = 0 .5. While in batch training a ﬁxed γ is suitable,
in a stochastic setting it is unlikely that a global decay par ameter is sufﬁcient. Our schedule has an
annealing effect in the sense that γ values near 1 are feasible late in training even with only 3-5 CG
iterations, a property that is otherwise hard to achieve. Th is allows us to beneﬁt from sharing more
information across iterations late in training, similar to that of a typical momentum method.
A remaining question to consider is how to set the sizes of the gradient and curvature mini-batches.
[24] discuss theoretical advantages to utilizing the same m ini-batches for computing the gradient
and curvature vector products. In our setting, this may lead to some difﬁculties. Using same-sized
batches allows λ →0 during training [24]. Unfortunately, this can become incom patible with our
short hard-limit on the number of CG iterations, since CG req uires more work to optimize ˆMθ
when λ approaches zero. To account for this, on classiﬁcation task s where 3-5 CG iterations are
used, we opt to use gradient mini-batches that are 5-10 times larger than curvature mini-batches.
For deep autoencoder tasks where more CG iterations are used , we instead set both gradient and
curvature batches to be the same size. The behavior of λ is dependent on whether or not dropout
is used during training. Figure 1 demonstrates the behavior of λ during classiﬁcation training with
and without the use of dropout. With dropout, λ no longer converges to 0 but instead plummets,
rises and ﬂattens out. In both settings, λ does not decrease substantially as to negatively effect the
proposed CG solution and consequently the reduction ratio. Thus, the amount of work required by
CG remains consistent late in training. The other beneﬁt to using larger gradient batches is to account
for the additional computation in computing curvature-vec tor products which would make training
longer if both mini-batches were small and of the same size. In [4], the gradients and objectives are
computed using the full training set throughout the algorit hm, including during CG backtracking
and the backtracking linesearch. We utilize the gradient mini-batch for the current iteration in order
to compute all necessary gradient and objectives throughout the algorithm.
4.2 Levenberg-Marquardt damping
Martens makes use of the following Levenberg-Marquardt style damping criteria for updating λ:
ifρ > 3
4, λ ←2
3λ elseifρ < 1
4, λ ←3
2 λ (13)
which given a suitable initial value will converge to zero as training progresses. We observed that
the above damping criteria is too harsh in the stochastic set ting in the sense that λ will frequently
oscillate, which is sensible given the size of the curvature mini-batches. We instead opt for a much
softer criterion, for which lambda is updated as
ifρ > 3
4, λ ← 99
100λ elseifρ < 1
4, λ ←100
99 λ (14)
This choice, although somewhat arbitrary, is consistently effective. Thus reduction ratio values
computed from curvature mini-batches will have less overall inﬂuence on the damping strength.
4.3 Integrating dropout
Dropout is a recently proposed method for improving the training of neural networks. During train-
ing, each hidden unit is omitted with a probability of 0.5 alo ng with optionally omitting input fea-
tures similar to that of a denoising autoencoder [25]. Dropo ut can be viewed in two ways. By
randomly omitting feature detectors, dropout prevents co-adaptation among detectors which can im-
prove generalization accuracy on held-out data. Secondly, dropout can be seen as a type of model
averaging. At test time, outgoing weights are halved. If we consider a network with a single hidden
layer and k feature detectors, using the mean network at test time corresponds to taking the geomet-
ric average of 2k networks with shared weights. Dropout is integrated in stochastic HF by randomly
omitting feature detectors on both gradient and curvature m ini-batches from the last hidden layer
5
0 100 200 300 400 5000
0.2
0.4
0.6
0.8
1
Epoch
lambda
 
 
0 200 400 600 800 10000
0.2
0.4
0.6
0.8
1
Epoch
lambda
 
 
Figure 1: Values of the damping strength λ during training of MNIST (left) and USPS (right) with
and without dropout using λ = 1 for classiﬁcation. When dropout is included, the damping strength
initially decreases followed by a steady increase over time.
during each iteration. Since we assume that the curvature mi ni-batches are a subset of the gradient
mini-batches, the same feature detectors are omitted in both cases.
Since the curvature estimates are noisy, it is important to c onsider the stability of updates when
different stochastic networks are used in each computation . The weight updates in dropout SGD
are augmented with momentum not only for stability but also t o speed up learning. Speciﬁcally, at
iteration k the parameter update is given by
∆ θk = pk∆ θk−1 −(1 −pk)αk⟨∇f⟩, θ k = θk−1 + ∆ θk (15)
where pk and ak are the momentum and learning rate, respectively. We incorp orate an additional
exponential decay term βe when performing parameter updates. Speciﬁcally, each parameter update
is computed as
θk = θk−1 + βeαkδk, β e = cβe−1 (16)
where c ∈(0, 1] is a ﬁxed parameter chosen by the user. Incorporating βe into the updates, along
with the use of δ-momentum, leads to more stable updates and ﬁne convergence particularly when
dropout is integrated during training.
4.4 Algorithm
Pseudo-code for one iteration of our implementation of stochastic Hessian-free is presented. Given
a gradient minibatch Xg
k and curvature minibatch Xc
k, we ﬁrst sample dropout units (if applicable)
for the inputs and last hidden layer of the network. These take the form of a binary vector, which are
multiplied component-wise by the activations yi. In our pseudo-code, CG (δ0
k, ∇f, P, ζ ) is used to
denote applying CG with initial solution δ0
k, gradient ∇f, pre-conditioner P and ζ iterations. Note
that, when computing δ-momentum, the ζ-th solution in iteration k −1 is used as opposed to the
solution chosen via backtracking. Given the objectivesfk−1 computed with θ and fk computed with
θ + δk, the reduction ratio ρ is calculated utilizing the un-damped quadratic approximation Mθ(δk).
This allows updating λ using the Levenberg-Marquardt style damping. Finally, a backtracking line-
search with at most ω steps is performed to compute the rate and serves as a last def ense against
potentially poor update directions.
Since curvature mini-batches are sampled from a subset of the gradient mini-batch, it is then sensible
to utilize different curvature mini-batches on different epochs. Along with cycling through gradient
mini-batches during each epoch, we also cycle through curvature subsets everyh epochs, where h is
the size of the gradient mini-batches divided by the size of the curvature mini-batches. For example,
if the gradient batch size is 1000 and the curvature batch siz e is 100, then curvature mini-batch
sampling completes a full cycle every 1000/100 = 10 epochs.
Finally, one simple way to speed up training as indicated in [ 24], is to cache the activations when
initially computing the objective fk. While each iteration of CG requires computing a curvature-
vector product, the network parameters are ﬁxed during CG an d is thus wasteful to re-compute the
network activations on each iteration.
6
Algorithm 1 Stochastic Hessian-Free Optimization
Xg
k ←gradient minibatch, Xc
k ←curvature minibatch, |Xg
k |= h|Xc
k|, h ∈Z+
Sample dropout units for inputs and last hidden layer
if start of new epoch then
γe ←min(1.01γe−1, .99) {δ-momentum}
end if
δ0
k ←γeδζ
k−1
fk−1 ←f(Xg
k ; θ), ∇f ←∇f(Xg
k ; θ), P ←Precon(Xg
k ; θ)
Solve (B + λI)δk = −∇f using CG(δ0
k, ∇f, P, ζ ) {Using Xc
k to compute Bδk}
fk ←f(Xg
k ; θ + δk) {CG backtracking}
for j = ζ - 1 to 1 do
f(θ + δj
k) ←f(Xg
k ; θ + δj
k)
if f(θ + δj
k) < f k then
fk ←f(θ + δj
k), δk ←δj
k
end if
end for
ρ ←(fk −fk−1)/(1
2 δT
k Bδk + ∇fT δk) {Using Xc
k to compute Bδk}
if ρ < . 25, λ ←1.01λ elseif ρ > . 75, λ ←.99λ end if
αk ←1, j ←0 {Backtracking linesearch}
while j < ω do
if fk > f k−1 + .01αk∇fT δk then αk ←.8αk, j ←j + 1 else break end if
end while
θ ←θ + βeαkδk, k ←k + 1 {Parameter update}
5 Experiments
We perform experimental evaluation on both classiﬁcation a nd deep autoencoder tasks. The goal
of classiﬁcation experiments is to determine the effective ness of SHF on test error generalization.
For autoencoder tasks, we instead focus just on measuring the effectiveness of the optimizer on the
training data. The datasets and experiments are summarized as follows:
•MNIST: Handwritten digits of size 28 ×28 with 60K training samples and 10K testing samples.
For classiﬁcation, we train networks of size 784-1200-1200-10 with rectiﬁer activations. For deep
autoencoders, the encoder architecture of 784-1000-500-250-30 with a symmetric decoding archi-
tecture is used. Logistic activations are used with a binary cross entropy error. For classiﬁcation
experiments, the data is scaled to have zero mean and unit variance.
•CURVES: Artiﬁcial dataset of curves of size 28 ×28 with 20K training samples and 10K testing
samples. We train a deep autoencoder using an encoding architecture of 784-400-200-100-50-25-
6 with symmetric decoding. Similar to MNIST, logistic activations and binary cross entropy error
are used.
•USPS: Handwritten digits of size 16 ×16 with 11K examples. We perform classiﬁcation using
5 randomly sampled batches of 8K training examples and 3K tes ting examples as in [26] Each
batch has an equal number of each digit. Classiﬁcation netwo rks of size 256-500-500-10 are
trained with rectiﬁer activations. The data is scaled to have zero mean and unit variance.
•Reuters: A collection of 8293 text documents from 65 categor ies. Each document is represented
as a 18900-dimensional bag-of-words vector. Word counts C are transformed to log( 1 + C)
as is done by [3]. The publically available train/test split of is used. We train networks of size
18900-65 for classiﬁcation due to the high dimensionality of the inputs, which reduces to softmax-
regression.
For classiﬁcation experiments, we perform comparison of SH F with and without dropout against
dropout SGD [3]. All classiﬁcation experiments utilize the sparse initialization of Martens [4] with
initial biases set to 0.1. The sparse initialization in comb ination with ReLUs make our networks
similar to the deep sparse rectiﬁer networks of [28]. All alg orithms are trained for 500 epochs on
MNIST and 1000 epochs on USPS and Reuters. We use weight decay of 5 ×10−4 for SHF and
2 ×10−5 for dropout SHF. A held-out validation set was used for determining the amount of input
7
0 100 200 300 400 5000
0.005
0.01
0.015
0.02
0.025
Epoch
classification error
MNIST
 
 
SHF
dSHF
dSGD−a
dSGD−l
0 200 400 600 800 10000
0.005
0.01
0.015
0.02
0.025
0.03
0.035
Epoch
classification error
USPS
 
 
SHF
dSHF
dSGD−a
dSGD−l
0 200 400 600 800 10000
0.02
0.04
0.06
0.08
0.1
0.12
0.14
0.16
0.18
0.2
Epoch
classification error
Reuters
 
 
SHF
dSHF
dSGD
Figure 2: Training and testing curves for classiﬁcation. dS HF: dropout SHF, dSGD: dropout SGD,
dSGD-a: dropout on all layers, dSGD-l: dropout on last hidde n layer only (as well as the inputs).
dropout for all algorithms. Both SHF and dropout SHF use initial damping of λ = 1 , gradient batch
size of 1000, curvature batch size of 100 and 3 CG iterations per batch.
Dropout SGD training uses an exponential decreasing learni ng rate schedule initialized at 10, in
combination with max-norm weight clipping [3]. This allows SGD to use larger learning rates for
greater exploration early in training. A linearly increasi ng momentum schedule is used with initial
momentum of 0.5 and ﬁnal momentum of 0.99. No weight decay is used. For additional comparison
we also train dropout SGD when dropout is only used in the last hidden layer, as is the case with
dropout SHF.
For deep autoencoder experiments, we use the same experimen tal setup as in Chapter 7 of [17].
In particular, we focus solely on training error without any L2 penalty in order to determine the
effectiveness of the optimizer on modeling the training dat a. Comparison is made against SGD,
SGD with momentum, HF and Nesterov’s accelerated gradient ( NAG). On CURVES, SHF uses an
initial damping of λ = 10 , gradient and curvature batch sizes of 2000 and 25 CG iterati ons per
batch. On MNIST, we use initial λ = 1 , gradient and curvature batch sizes of 3000 and 50 CG
iterations per batch. Autoencoder training is ran until no s ufﬁcient progress is made, which occurs
at around 250 epochs on CURVES and 100 epochs on MNIST.
5.1 Classiﬁcation results
Figure 2 summarizes our classiﬁcation results. At epoch 500, dropou t SHF achieves 107 errors on
MNIST. This result is similar to [3] which achieve 100-115 errors with various network sizes when
training for a few thousand epochs. Without dropout or input corruption, SHF achieves 159 errors
on MNIST, on par with existing methods that do not incorporat e prior knowledge, pre-training,
image distortions or dropout. As with [4], we hypothesize th at further improvements can be made
by ﬁne-tuning with SHF after unsupervised layerwise pre-training.
After 1000 epochs of training on ﬁve random splits of USPS, we obtain ﬁnal classiﬁcation errors of
1%, 1.1%, 0.8%, 0.9% and 0.97% with a mean test error of 0.95%. Both algorithms use 50% input
corruption. For additional comparison, [29] obtains a mean classiﬁcation error of 1.14% using a
pre-trained deep network for large-margin nearest neighbo r classiﬁcation with the same size splits.
Without dropout, SHF overﬁts the training data.
On the Reuters dataset, SHF with and without dropout both dem onstrate accelerated training. We
hypothesize that further speedup may also be obtained by sta rting training with a much smaller λ
initialization, which we suspect is conservative given that the problem is convex.
8
Table 1: Training errors on the deep autoencoder tasks. All results are obtained from [17]. M(0.99)
refers to momentum capped at 0.99 and similarily for M(0.9). SGD-VI refers to SGD using the
variance normalized initialization of [15].
problem NAG M(0.99) M(0.9) SGD SGD-VI [19] HF SHF
CURVES 0.078 0.110 0.220 0.250 0.160 0.110 0.089
MNIST 0.730 0.770 0.990 1.100 0.900 0.780 0.877
0 50 100 150 200 2500
0.05
0.1
0.15
0.2
0.25
Epoch
train_L2
CURVES
 
 
SHF
NAG
HF
SGD−VI
SGD
0 20 40 60 80 1000
0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.8
Epoch
train_L2
MNIST
 
 
SHF
NAG
HF
SGD−VI
SGD
Figure 3: Learning curves for the deep autoencoder tasks. Th e CG decay parameter γ is shut off at
epoch 160 on CURVES and epoch 60 on MNIST.
5.2 Deep autoencoder results
Figure 3 and table 1 summarize our results. Inspired by [17] we make one addition al modiﬁcation
to our algorithms. As soon as training begins to diverge, we t urn off the CG decay parameter γ
in a similar fashion as the the momentum parameter µ is decreased in [17]. When γ = 0 , CG
is no longer initialized from the previous solution and is in stead initialized to zero. As with [17],
this has a dramatic effect on the training error but to a lesse r extent as momentum and Nesterov’s
accelerated gradient. [17] describes the behaviour of this effect as follows: with a large momentum,
the optimizer is able to make steady progress along slow chan ging directions of low curvature.
By decreasing the momentum late in training, the optimizer i s then able to quickly reach a local
minimum from ﬁner optimization along high curvature direct ions, which would otherwise be too
difﬁcult to obtain with an aggressive momentum schedule. Th is observation further motivates the
relationship between momentum and information sharing through CG.
Our experimental results demonstrate that SHF does not perf orm signiﬁcantly better or worse on
these datasets compared to existing approaches. It is able t o outperform HF on CURVES but not
on MNIST. An attractive property that is shared with both HF a nd SHF is not requiring the careful
schedule tuning that is necessary for momentum and NAG. We also attempted experiments with SHF
using the same setup for classiﬁcation with smaller batches and 5 CG iterations. The results were
worse: on CURVES the lowest training error obtained was 0.19. This shows that while such a setup
is useful from the viewpoint of noisy updates and test generalization, they hamper the effectiveness
of making progress on hard to optimize regions.
6 Conclusion
In this paper we proposed a stochastic variation of Martens’Hessian-free optimization incorporating
dropout for training neural networks on classiﬁcation and d eep autoencoder tasks. By adapting the
batch sizes and number of CG iterations, SHF can be construct ed to perform well for classiﬁcation
9
against dropout SGD or optimizing deep autoencoders compar ing HF, NAG and momentum meth-
ods. While our initial results are promising, of interest wo uld be adapting stochastic Hessian-free
optimization to other network architectures:
•Convolutional networks. The most common approach to training convolutional network s has
been SGD incorporating a diagonal Hessian approximation [8]. Dropout SGD was recently used
for training a deep convolutional network on ImageNet [30].
•Recurrent Networks. It was largely believed that RNNs were too difﬁcult to train with SGD due
to the exploding/vanishing gradient problem. In recent yea rs, recurrent networks have become
popular again due to several advancements made in their training [31].
•Recursive Networks. Recursive networks have been successfully used for tasks such as sentiment
classiﬁcation and compositional modeling of natural language from word embeddings [32]. These
architectures are usually trained using L-BFGS.
It is not clear yet whether this setup is easily generalizabl e to the above architectures or whether
improvements need to be considered. Furthermore, addition al experimental comparison would in-
volve dropout SGD with the adaptive methods of Adagrad [9] or [11], as well as the importance of
pre-conditioning CG. None the less, we hope that this work in itiates future research in developing
stochastic Hessian-free algorithms.
Acknowledgments
The author would like to thank Csaba Szepesvári for helpful discussion as well as David Sussillo for
his guidance when ﬁrst learning about and implementing HF. T he author would also like to thank
the anonymous ICLR reviewers for their comments and suggestions.
References
[1] L. Bottou and O. Bousquet. The tradeoffs of large-scale l earning. Optimization for Machine
Learning, page 351, 2011.
[2] Y . Bengio, P. Lamblin, D. Popovici, and H. Larochelle. Gr eedy layer-wise training of deep
networks. NIPS, 19:153, 2007.
[3] G.E. Hinton, N. Srivastava, A. Krizhevsky, I. Sutskever , and R.R. Salakhutdinov. Improving
neural networks by preventing co-adaptation of feature detectors. arXiv:1207.0580, 2012.
[4] J. Martens. Deep learning via hessian-free optimizatio n. In ICML, volume 951, 2010.
[5] J. Martens and I. Sutskever. Learning recurrent neural networks with hessian-free optimization.
In ICML, 2011.
[6] B.A. Pearlmutter. Fast exact multiplication by the hessian. Neural Computation, 6(1):147–160,
1994.
[7] N.N. Schraudolph. Fast curvature matrix-vector produc ts for second-order gradient descent.
Neural computation, 14(7):1723–1738, 2002.
[8] Y . LeCun, L. Bottou, G. Orr, and K. Müller. Efﬁcient backp rop. Neural networks: T ricks of
the trade, pages 546–546, 1998.
[9] J. Duchi, E. Hazan, and Y . Singer. Adaptive subgradient m ethods for online learning and
stochastic optimization. JMLR, 12:2121–2159, 2010.
[10] J. Dean, G. Corrado, R. Monga, K. Chen, M. Devin, Q. Le, M. Mao, A. Senior, P. Tucker,
K. Yang, et al. Large scale distributed deep networks. In NIPS, pages 1232–1240, 2012.
[11] T. Schaul, S. Zhang, and Y . LeCun. No more pesky learning rates. arXiv:1206.1106, 2012.
[12] A. Bordes, L. Bottou, and P. Gallinari. Sgd-qn: Careful quasi-newton stochastic gradient
descent. JMLR, 10:1737–1754.
[13] Razvan Pascanu and Yoshua Bengio. Natural gradient rev isited. arXiv preprint
arXiv:1301.3584, 2013.
[14] N. Le Roux, P.A. Manzagol, and Y . Bengio. Topmoumoute online natural gradient algorithm.
In NIPS, 2007.
10
[15] Xavier Glorot and Yoshua Bengio. Understanding the dif ﬁculty of training deep feedforward
neural networks. In AISTATS, 2010.
[16] Yann N Dauphin and Yoshua Bengio. Big neural networks wa ste capacity. arXiv preprint
arXiv:1301.3583, 2013.
[17] I. Sutskever. T raining Recurrent Neural Networks. PhD thesis, University of Toronto, 2013.
[18] R.H. Byrd, G.M. Chin, J. Nocedal, and Y . Wu. Sample size selection in optimization methods
for machine learning. Mathematical Programming, pages 1–29, 2012.
[19] O. Chapelle and D. Erhan. Improved preconditioner for h essian free optimization. In NIPS
W orkshop on Deep Learning and Unsupervised F eature Learning, 2011.
[20] I. Sutskever, J. Martens, and G. Hinton. Generating tex t with recurrent neural networks. In
ICML, 2011.
[21] O. Vinyals and D. Povey. Krylov subspace descent for dee p learning. arXiv:1111.4259, 2011.
[22] Q.V . Le, J. Ngiam, A. Coates, A. Lahiri, B. Prochnow, andA.Y . Ng. On optimization methods
for deep learning. In ICML, 2011.
[23] Y . Bengio. Practical recommendations for gradient-ba sed training of deep architectures.
arXiv:1206.5533, 2012.
[24] J. Martens and I. Sutskever. Training deep and recurren t networks with hessian-free optimiza-
tion. Neural Networks: T ricks of the T rade, pages 479–535, 2012.
[25] P. Vincent, H. Larochelle, Y . Bengio, and P.A. Manzagol . Extracting and composing robust
features with denoising autoencoders. In ICML, pages 1096–1103, 2008.
[26] R. Min, D.A. Stanley, Z. Yuan, A. Bonner, and Z. Zhang. A d eep non-linear feature mapping
for large-margin knn classiﬁcation. In Ninth IEEE International Conference on Data Mining ,
pages 357–366. IEEE, 2009.
[27] Deng Cai, Xuanhui Wang, and Xiaofei He. Probabilistic d yadic data analysis with local and
global consistency. In ICML, pages 105–112, 2009.
[28] X. Glorot, A. Bordes, and Y . Bengio. Deep sparse rectiﬁer neural networks. In AISTATS, 2011.
[29] R. Min, D.A. Stanley, Z. Yuan, A. Bonner, and Z. Zhang. A d eep non-linear feature mapping
for large-margin knn classiﬁcation. In ICDM, pages 357–366, 2009.
[30] A. Krizhevsky, I. Sutskever, and G. Hinton. Imagenet cl assiﬁcation with deep convolutional
neural networks. NIPS, 25, 2012.
[31] Y . Bengio, N. Boulanger-Lewandowski, and R. Pascanu. A dvances in optimizing recurrent
networks. arXiv:1212.0901, 2012.
[32] R. Socher, B. Huval, C.D. Manning, and A.Y . Ng. Semanticcompositionality through recursive
matrix-vector spaces. In EMNLP, pages 1201–1211, 2012.
11


Peer review:
[Review]:
I want to say thanks again to the conference organizers, reviewers and openreview.net developers for doing a great job.

I have updated the code on my webpage to include two additional features: max norm weight clipping and training deep autoencoders. Autoencoder training uses symmetric encoding / decoding and supports denoising and L2 penalties.