Training Quantized Nets: A Deeper Understanding

arXiv:1706.02379v3 [cs.LG] 13 Nov 2017

Hao Li1∗, Soham De1∗, Zheng Xu1 , Christoph Studer2 , Hanan Samet1 , Tom Goldstein1
1
Department of Computer Science, University of Maryland, College Park
2
School of Electrical and Computer Engineering, Cornell University
{haoli,sohamde,xuzh,hjs,tomg}@cs.umd.edu, studer@cornell.edu

Abstract
Currently, deep neural networks are deployed on low-power portable devices by first training
a full-precision model using powerful hardware, and then deriving a corresponding lowprecision model for efficient inference on such systems. However, training models directly
with coarsely quantized weights is a key step towards learning on embedded platforms that
have limited computing resources, memory capacity, and power consumption. Numerous
recent publications have studied methods for training quantized networks, but these studies
have mostly been empirical. In this work, we investigate training methods for quantized neural networks from a theoretical viewpoint. We first explore accuracy guarantees for training
methods under convexity assumptions. We then look at the behavior of these algorithms for
non-convex problems, and show that training algorithms that exploit high-precision representations have an important greedy search phase that purely quantized training methods
lack, which explains the difficulty of training using low-precision arithmetic.

1

Introduction

Deep neural networks are an integral part of state-of-the-art computer vision and natural language
processing systems. Because of their high memory requirements and computational complexity,
networks are usually trained using powerful hardware. There is an increasing interest in training
and deploying neural networks directly on battery-powered devices, such as cell phones or other
platforms. Such low-power embedded systems are memory and power limited, and in some cases
lack basic support for floating-point arithmetic.
To make neural nets practical on embedded systems, many researchers have focused on training nets
with coarsely quantized weights. For example, weights may be constrained to take on integer/binary
values, or may be represented using low-precision (8 bits or less) fixed-point numbers. Quantized nets
offer the potential of superior memory and computation efficiency, while achieving performance that
is competitive with state-of-the-art high-precision nets. Quantized weights can dramatically reduce
memory size and access bandwidth, increase power efficiency, exploit hardware-friendly bitwise
operations, and accelerate inference throughput [1–3].
Handling low-precision weights is difficult and motivates interest in new training methods. When
learning rates are small, stochastic gradient methods make small updates to weight parameters.
Binarization/discretization of weights after each training iteration “rounds off” these small updates
and causes training to stagnate [1]. Thus, the naïve approach of quantizing weights using a rounding
procedure yields poor results when weights are represented using a small number of bits. Other
approaches include classical stochastic rounding methods [4], as well as schemes that combine
full-precision floating-point weights with discrete rounding procedures [5]. While some of these
schemes seem to work in practice, results in this area are largely experimental, and little work has
been devoted to explaining the excellent performance of some methods, the poor performance of
others, and the important differences in behavior between these methods.
∗

Equal contribution. Author ordering determined by a cryptographically secure random number generator.

31st Conference on Neural Information Processing Systems (NIPS 2017), Long Beach, CA, USA.

Contributions This paper studies quantized training methods from a theoretical perspective, with
the goal of understanding the differences in behavior, and reasons for success or failure, of various
methods. In particular, we present a convergence analysis showing that classical stochastic rounding
(SR) methods [4] as well as newer and more powerful methods like BinaryConnect (BC) [5] are
capable of solving convex discrete problems up to a level of accuracy that depends on the quantization
level. We then address the issue of why algorithms that maintain floating-point representations, like
BC, work so well, while fully quantized training methods like SR stall before training is complete.
We show that the long-term behavior of BC has an important annealing property that is needed for
non-convex optimization, while classical rounding methods lack this property.

2

Background and Related Work

The arithmetic operations of deep networks can be truncated down to 8-bit fixed-point without
significant deterioration in inference performance [4, 6–9]. The most extreme scenario of quantization
is binarization, in which only 1-bit (two states) is used for weight representation [10, 5, 1, 3, 11, 12].
Previous work on obtaining a quantized neural network can be divided into two categories: quantizing
pre-trained models with or without retraining [7, 13, 6, 14, 15], and training a quantized model from
scratch [4, 5, 3, 1, 16]. We focus on approaches that belong to the second category, as they can be
used for both training and inference under constrained resources.
For training quantized NNs from scratch, many authors suggest maintaining a high-precision floating
point copy of the weights while feeding quantized weights into backprop [5, 11, 3, 16], which results
in good empirical performance. There are limitations in using such methods on low-power devices,
however, where floating-point arithmetic is not always available or not desirable. Another widely
used solution using only low-precision weights is stochastic rounding [17, 4]. Experiments show
that networks using 16-bit fixed-point representations with stochastic rounding can deliver results
nearly identical to 32-bit floating-point computations [4], while lowering the precision down to 3-bit
fixed-point often results in a significant performance degradation [18]. Bayesian learning has also
been applied to train binary networks [19, 20]. A more comprehensive review can be found in [3].

3

Training Quantized Neural Nets

We consider empirical risk minimization problems of the form:
m
1 X
min F (w) :=
fi (w),
w∈W
m i=1

(1)

where the objective function decomposes into a sum over many functions fi : Rd → R. Neural
networks have objective functions of this form where each fi is a non-convex loss function. When
floating-point representations are available, the standard method for training neural networks is
stochastic gradient descent (SGD), which on each iteration selects a function f˜ randomly from
{f1 , f2 , . . . , fm }, and then computes
SGD: wt+1 = wt − αt ∇f˜(wt ),
(2)
for some learning rate αt . In this paper, we consider the problem of training convolutional neural
networks (CNNs). Convolutions are computationally expensive; low precision weights can be used
to accelerate them by replacing expensive multiplications with efficient addition and subtraction
operations [3, 9] or bitwise operations [11, 16].
To train networks using a low-precision representation of the weights, a quantization function Q(·)
is needed to convert a real-valued number w into a quantized/rounded version ŵ = Q(w). We use
the same notation for quantizing vectors, where we assume Q acts on each dimension of the vector.
Different quantized optimization routines can be defined by selecting different quantizers, and also
by selecting when quantization happens during optimization. The common options are:
Deterministic Rounding (R) A basic uniform or deterministic quantization function snaps a
floating point value to the closest quantized value as:


|w| 1
Qd (w) = sign(w) · ∆ ·
+
,
(3)
∆
2
2

where ∆ denotes the quantization step or resolution, i.e., the smallest positive number that is
representable. One exception to this definition is when we consider binary weights, where all weights
are constrained to have two values w ∈ {−1, 1} and uniform rounding becomes Qd (w) = sign(w).
The deterministic rounding SGD maintains quantized weights with updates of the form:

Deterministic Rounding: wbt+1 = Qd wbt − αt ∇f˜(wbt ) ,

(4)

where wb denotes the low-precision weights, which are quantized using Qd immediately after applying
the gradient descent update. If gradient updates are significantly smaller than the quantization step,
this method loses gradient information and weights may never be modified from their starting values.
Stochastic Rounding (SR) The quantization function for stochastic rounding is defined as:
 w
w
w
b ∆ c + 1 for p ≤ ∆
− b∆
c,
Qs (w) = ∆ ·
w
b∆c
otherwise,

(5)

where p ∈ [0, 1] is produced by a uniform random number generator. This operator is nondeterministic, and rounds its argument up with probability w/∆ − bw/∆c, and down otherwise.
This quantizer satisfies the important property E[Qs (w)] = w. Similar to the deterministic rounding
method, the SR optimization method also maintains quantized weights with updates of the form:

Stochastic Rounding: wbt+1 = Qs wbt − αt ∇f˜(wbt ) .
(6)
BinaryConnect (BC) The BinaryConnect algorithm [5] accumulates gradient updates using a
full-precision buffer wr , and quantizes weights just before gradient computations as follows.

BinaryConnect: wrt+1 = wrt − αt ∇f˜ Q(wrt ) .
(7)
Either stochastic rounding Qs or deterministic rounding Qd can be used for quantizing the weights
wr , but in practice, Qd is the common choice. The original BinaryConnect paper constrains the
low-precision weights to be {−1, 1}, which can be generalized to {−∆, ∆}. A more recent method,
Binary-Weights-Net (BWN) [3], allows different filters to have different scales for quantization,
which often results in better performance on large datasets.
Notation For the rest of the paper, we use Q to denote both Qs and Qd unless the situation requires
this to be distinguished. We also drop the subscripts on wr and wb , and simply write w.

4

Convergence Analysis

We now present convergence guarantees for the Stochastic Rounding (SR) and BinaryConnect
(BC) algorithms, with updates of the form (6) and (7), respectively. For the purposes of deriving
theoretical guarantees, we assume each fi in (1) is differentiable and convex, and the domain
W is convex and has dimension d. We consider both the case where F is µ-strongly convex:
h∇F (w0 ), w − w0 i ≤ F (w) − F (w0 ) − µ2 kw − w0 k2 , as well as where F is weakly convex. We also
assume the (stochastic) gradients are bounded: Ek∇f˜(wt )k2 ≤ G2 . Some results below also assume
the domain of the problem is finite. In this case, the rounding algorithm clips values that leave the
domain. For example, in the binary case, rounding returns bounded values in {−1, 1}.
4.1

Convergence of Stochastic Rounding (SR)

We can rewrite the update rule (6) as:
wt+1 = wt − αt ∇f˜(wt ) + rt ,
where rt = Qs (wt − αt ∇f˜(wt )) − wt + αt ∇f˜(wt ) denotes the quantization error on the t-th
iteration. We want to bound this error in expectation. To this end, we present the following lemma.
Lemma 1. The stochastic rounding error rt on each iteration can be bounded, in expectation, as:
√
2
E rt ≤ d∆αt G,
where d denotes the dimension of w.
3

Proofs for all theoretical results are presented in the Appendices. From Lemma 1, we see that
the rounding error per step decreases as the learning rate αt decreases. This is intuitive since the
probability of an entry in wt+1 differing from wt is small when the gradient update is small relative
to ∆. Using the above lemma, we now present convergence rate results for Stochastic Rounding (SR)
in both the strongly-convex case and the non-strongly convex case. Our error estimates are ergodic,
PT
i.e., they are in terms of w̄T = T1 t=1 wt , the average of the iterates.
Theorem 1. Assume that F is µ-strongly convex and the learning rates are given by αt =
Consider the SR algorithm with updates of the form (6). Then, we have:
√
(1 + log(T + 1))G2
d∆G
T
?
E[F (w̄ ) − F (w )] ≤
+
,
2µT
2
where w? = arg minw F (w).
Theorem 2. Assume the domain has finite diameter D, and learning rates are given by αt =
a constant c. Consider the SR algorithm with updates of the form (6). Then, we have:
√
√
1
T +1 2
d∆G
T
?
2
cG +
.
E[F (w̄ ) − F (w )] ≤ √ D +
2T
2
c T

1
µ(t+1) .

c
√
,
t

for

We see that in both cases, SR converges until it reaches an “accuracy floor.” As the quantization
becomes more fine grained, our theory predicts that the accuracy of SR approaches that of highprecision floating point at a rate linear in ∆. This extra term caused by the discretization is unavoidable
since this method maintains quantized weights.
4.2

Convergence of Binary Connect (BC)

When analyzing the BC algorithm, we assume that the Hessian satisfies the Lipschitz bound:
k∇2 fi (x) − ∇2 fi (y)k ≤ L2 kx − yk for some L2 ≥ 0. While this is a slightly non-standard
assumption, we will see that it enables us to gain better insights into the behavior of the algorithm.
The results here hold for both stochastic and uniform rounding. In this case, the quantization error r
does not approach 0 as in SR-SGD. Nonetheless, the effect of this rounding error diminishes with
shrinking αt because αt multiplies the gradient update, and thus implicitly the rounding error as well.
Theorem 3. Assume F is L-Lipschitz smooth, the domain has finite diameter D, and learning rates
are given by αt = √ct . Consider the BC-SGD algorithm with updates of the form (7). Then, we have:
√
T +1 2 √
1
T
?
2
E[F (w̄ ) − F (w )] ≤ √ D +
cG + d∆LD.
2T
2c T
As with SR, BC can only converge up to an error floor. So far this looks a lot like the convergence
guarantees for SR. However, things change when we assume strong convexity and bounded Hessian.
1
Theorem 4. Assume that F is µ-strongly convex and the learning rates are given by αt = µ(t+1)
.
Consider the BC algorithm with updates of the form (7). Then we have:
√
(1 + log(T + 1))G2
DL2 d∆
+
.
E[F (w̄T ) − F (w? )] ≤
2µT
2
Now, the error floor is determined by both ∆ and L2 . For a quadratic least-squares problem, the
gradient of F is linear and the Hessian is constant. Thus, L2 = 0 and we get the following corollary.
1
Corollary 1. Assume that F is quadratic and the learning rates are given by αt = µ(t+1)
. The BC
algorithm with updates of the form (7) yields
E[F (w̄T ) − F (w? )] ≤

(1 + log(T + 1))G2
.
2µT

We see that the real-valued weights accumulated in BC can converge to the true minimizer of quadratic
losses. Furthermore, this suggests that, when the function behaves like a quadratic on the distance
4

Figure 1: The SR method starts at some location x (in this case 0), adds a perturbation to x, and then rounds.
As the learning rate α gets smaller, the distribution of the perturbation gets “squished” near the origin, making
the algorithm less likely to move. The “squishing” effect is the same for the part of the distribution lying to the
left and to the right of x, and so it does not effect the relative probability of moving left or right.
scale ∆, one would expect BC to perform fundamentally better than SR. While this may seem
like a restrictive condition, there is evidence that even non-convex neural networks become well
approximated as a quadratic in the later stages of optimization within a neighborhood of a local
minimum [21].
Note, our convergence results on BC are for wr instead of wb , and these measures of convergence are
not directly comparable. It is not possible to bound wb when BC is used, as the values of wb may
not converge in the usual sense (e.g., in the +/-1 binary case wr might converge to 0, in which case
arbitrarily small perturbations to wr might send wb to +1 or -1).

5

What About Non-Convex Problems?

The global convergence results presented above for convex problems show that, in general, both
the SR and BC algorithms converge to within O(∆) accuracy of the minimizer (in expected value).
However, these results do not explain the large differences between these methods when applied to
non-convex neural nets. We now study how the long-term behavior of SR differs from BC. Note
that this section makes no convexity assumptions, and the proposed theoretical results are directly
applicable to neural networks.
Typical (continuous-valued) SGD methods have an important exploration-exploitation tradeoff. When
the learning rate is large, the algorithm explores by moving quickly between states. Exploitation
happens when the learning rate is small. In this case, noise averaging causes the algorithm more
greedily pursues local minimizers with lower loss values. Thus, the distribution of iterates produced
by the algorithm becomes increasingly concentrated near minimizers as the learning rate vanishes
(see, e.g., the large-deviation estimates in [22]). BC maintains this property as well—indeed, we saw
in Corollary 1 a class of problems for which the iterates concentrate on the minimizer for small αt .
In this section, we show that the SR method lacks this important tradeoff: as the stepsize gets small
and the algorithm slows down, the quality of the iterates produced by the algorithm does not improve,
and the algorithm does not become progressively more likely to produce low-loss iterates. This
behavior is illustrated in Figures 1 and 2.
To understand this problem conceptually, consider the simple case of a one-variable optimization
problem starting at x0 = 0 with ∆ = 1 (Figure 1). On each iteration, the algorithm computes a
stochastic approximation ∇f˜ of the gradient by sampling from a distribution, which we call p. This
gradient is then multiplied by the stepsize to get α∇f˜. The probability of moving to the right (or
left) is then roughly proportional to the magnitude of α∇f˜. Note the random variable α∇f˜ has
distribution pα (z) = α−1 p(z/α).
Now, suppose that α is small enough that we can neglect the tails of pα (z) that lie outside the interval
[−1, 1]. The probability of transitioning from x0 = 0 to x1 = 1 using stochastic rounding, denoted
by Tα (0, 1), is then
Z 1
Z
Z 1/α
Z ∞
1 1
Tα (0, 1) ≈
zpα (z)dz =
zp(z/α) dz = α
p(x)x dx ≈ α
p(x)x dx,
α 0
0
0
0
where the first approximation is because we neglected the unlikely case that α∇f˜ > 1, and the
second approximation appears because we added a small tail probability to the estimate. These
5

12

Loss Value

10
8
6
4
2
0

-2

0

2

4

6

8

Weight w

(a) α = 1.0

(b) α = 0.1

(c) α = 0.01

(d) α = 0.001

Figure 2: Effect of shrinking the learning rate in SR vs BC on a toy problem. The left figure plots the objective

function (8). Histograms plot the distribution of the quantized weights over 106 iterations. The top row of plots
correspond to BC, while the bottom row is SR, for different learning rates α. As the learning rate α shrinks, the
BC distribution concentrates on a minimizer, while the SR distribution stagnates.

approximations get more accurate for small α. We see that, assuming the tails of p are “light” enough,
R∞
R0
we have Tα (0, 1) ∼ α 0 p(x)x dx as α → 0. Similarly, Tα (0, −1) ∼ α −∞ p(x)x dx as α → 0.
What does this observation mean for the behavior of SR? First of all, the probability of leaving x0 on
an iteration is
Z ∞

Z 0
Tα (0, −1) + Tα (0, 1) ≈ α
p(x)x dx +
p(x)x dx ,
−∞

0

which vanishes for small α. This means the algorithm slows down as the learning rate drops off,
which is not surprising. However, the conditional probability of ending up at x1 = 1 given that the
algorithm did leave x0 is
R∞
p(x)x dx
Tα (0, 1)
1
0
0
= R0
Tα (0, 1|x 6= x ) ≈
,
R∞
Tα (0, −1) + Tα (0, 1)
p(x)x dx +
p(x)x dx
−∞

0

which does not depend on α. In other words, provided α is small, SR, on average, makes the same
decisions/transitions with learning rate α as it does with learning rate α/10; it just takes 10 times
longer to make those decisions when α/10 is used. In this situation, there is no exploitation benefit in
decreasing α.
5.1

Toy Problem

To gain more intuition about the effect of shrinking the learning rate in SR vs BC, consider the
following simple 1-dimensional non-convex problem:
 2
if w < 1,
w + 2,
min f (w) := (w − 2.5)2 + 0.75,
(8)
if 1 ≤ w < 3.5,
w

(w − 4.75)2 + 0.19, if w ≥ 3.5.
Figure 2 shows a plot of this loss function. To visualize the distribution of iterates, we initialize at
w = 4.0, and run SR and BC for 106 iterations using a quantization resolution of 0.5.
Figure 2 shows the distribution of the quantized weight parameters w over the iterations when
optimized with SR and BC for different learning rates α. As we shift from α = 1 to α = 0.001, the
distribution of BC iterates transitions from a wide/explorative distribution to a narrow distribution
in which iterates aggressively concentrate on the minimizer. In contrast, the distribution produced
by SR concentrates only slightly and then stagnates; the iterates are spread widely even when the
learning rate is small.
5.2

Asymptotic Analysis of Stochastic Rounding

The above argument is intuitive, but also informal. To make these statements rigorous, we interpret
the SR method as a Markov chain. On each iteration, SR starts at some state (iterate) x, and moves to
6

0.6

0.2

0.8

0.6

0.4

0.2
B

A

B

A

0.2
0.6

0.1
0.3

0.4
0.2

0.2

0.2
0.1

0.1

C

C

0.2

0.6

Figure 3: Markov chain example with 3 states. In the right figure, we halved each transition probability for
moving between states, with the remaining probability put on the self-loop. Notice that halving all the transition
probabilities would not change the equilibrium distribution, and instead would only increase the mixing time of
the Markov chain.

a new state y with some transition probability Tα (x, y) that depends only on x and the learning rate
α. For fixed α, this is clearly a Markov process with transition matrix2 Tα (x, y).
The long-term behavior of this Markov process is determined by the stationary distribution of
Tα (x, y). We show below that for small α, the stationary distribution of Tα (x, y) is nearly invariant
to α, and thus decreasing α below some threshold has virtually no effect on the long term behavior of
the method. This happens because, as α shrinks, the relative transition probabilities remain the same
(conditioned on the fact that the parameters change), even though the absolute probabilities decrease
(see Figure 3). In this case, there is no exploitation benefit to decreasing α.
Theorem 5. Let px,k denote the probability distribution of the kth entry in ∇f˜(x), the stochastic gradient estimate at x. Assume there is a constant C1 such that for all x, k, and ν we have
R∞
R C2
R0
1
px,k (z) dz ≤ C
ν 2 , and some C2 such that both 0 px,k (z) dz > 0 and −C2 px,k (z) dz > 0.
ν
Define the matrix
R ∞
z

R0 px,k (z) ∆ dz, if x and y differ only at coordinate k, and yk = xk + ∆
0
Ũ (x, y) =
p (z) z dz, if x and y differ only at coordinate k, and yk = xk − ∆
 −∞ x,k ∆

0, otherwise,
and the associated markov chain transition matrix
T̃α0 = I − α0 · diag(1T Ũ ) + α0 Ũ ,

(9)

where α0 is the largest constant that makes T̃α0 non-negative. Suppose T̃α has a stationary distribution, denoted π̃. Then, for sufficiently small α, Tα has a stationary distribution πα , and
lim πα = π̃.

α→0

Furthermore, this limiting distribution satisfies π̃(x) > 0 for any state x, and is thus not concentrated
on local minimizers of f .
While the long term stationary behavior of SR is relatively insensitive to α, the convergence speed
of the algorithm is not. To measure this, we consider the mixing time of the Markov chain. Let πα
denote the stationary distribution of a Markov chain. We say that the -mixing time of the chain is
M if M is the smallest integer such that [23]
|P(xM ∈ A|x0 ) − π(A)| ≤ , for all x0 and all subsets of states A ⊆ X.

(10)

We show below that the mixing time of the Markov chain gets large for small α, which means
exploration slows down, even though no exploitation gain is being realized.
Theorem 6. Let px,k satisfy the assumptions of Theorem 5. Choose some  sufficiently small that
there exists a proper subset of states A ⊂ X with stationary probability πα (A) greater than . Let
M (α) denote the -mixing time of the chain with learning rate α. Then,
lim M (α) = ∞.

α→0
2

Our analysis below does not require the state space to be finite, so Tα (x, y) may be a linear operator rather
than a matrix. Nonetheless, we use the term “matrix” as it is standard.

7

Table 1: Top-1 test error after training with full-precision (ADAM), binarized weights (R-ADAM, SR-ADAM,
BC-ADAM), and binarized weights with big batch size (Big SR-ADAM).

CIFAR-10

CIFAR-100 ImageNet

VGG-9 VGG-BC ResNet-56 WRN-56-2 ResNet-56 ResNet-18
ADAM
BC-ADAM
Big SR-ADAM
SR-ADAM
R-ADAM

6

7.97
10.36
16.95
23.33
23.99

7.12
8.21
16.77
20.56
21.88

8.10
8.83
19.84
26.49
33.56

6.62
7.17
16.04
21.58
27.90

33.98
35.34
50.79
58.06
68.39

36.04
52.11
77.68
88.86
91.07

Experiments

To explore the implications of the theory above, we train both VGG-like networks [24] and Residual
networks [25] with binarized weights on image classification problems. On CIFAR-10, we train
ResNet-56, wide ResNet-56 (WRN-56-2, with 2X more filters than ResNet-56), VGG-9, and the
high capacity VGG-BC network used for the original BC model [5]. We also train ResNet-56 on
CIFAR-100, and ResNet-18 on ImageNet [26].
We use Adam [27] as our baseline optimizer as we found it to frequently give better results than
well-tuned SGD (an observation that is consistent with previous papers on quantized models [1–5]),
and we train with the three quantized algorithms mentioned in Section 3, i.e., R-ADAM, SR-ADAM
and BC-ADAM. The image pre-processing and data augmentation procedures are the same as [25].
Following [3], we only quantize the weights in the convolutional layers, but not linear layers, during
training (See Appendix H.1 for a discussion of this issue, and a detailed description of experiments).
We set the initial learning rate to 0.01 and decrease the learning rate by a factor of 10 at epochs 82 and
122 for CIFAR-10 and CIFAR-100 [25]. For ImageNet experiments, we train the model for 90 epochs
and decrease the learning rate at epochs 30 and 60. See Appendix H for additional experiments.
Results The overall results are summarized in Table 1. The binary model trained by BC-ADAM
has comparable performance to the full-precision model trained by ADAM. SR-ADAM outperforms
R-ADAM, which verifies the effectiveness of Stochastic Rounding. There is a performance gap
between SR-ADAM and BC-ADAM across all models and datasets. This is consistent with our
theoretical results in Sections 4 and 5, which predict that keeping track of the real-valued weights as
in BC-ADAM should produce better minimizers.
Exploration vs exploitation tradeoffs Section 5 discusses the exploration/exploitation tradeoff
of continuous-valued SGD methods and predicts that fully discrete methods like SR are unable to
enter a greedy phase. To test this effect, we plot the percentage of changed weights (signs different
from the initialization) as a function of the training epochs (Figures 4 and 5). SR-ADAM explores
aggressively; it changes more weights in the conv layers than both R-ADAM and BC-ADAM, and
keeps changing weights until nearly 40% of the weights differ from their starting values (in a binary
model, randomly re-assigning weights would result in 50% change). The BC method never changes
more than 20% of the weights (Fig 4(b)), indicating that it stays near a local minimizer and explores
less. Interestingly, we see that the weights of the conv layers were not changed at all by R-ADAM;
when the tails of the stochastic gradient distribution are light, this method is ineffective.
6.1

A Way Forward: Big Batch Training

We saw in Section 5 that SR is unable to exploit local minima because, for small learning rates,
shrinking the learning rate does not produce additional bias towards moving downhill. This was
illustrated in Figure 1. If this is truly the cause of the problem, then our theory predicts that we can
improve the performance of SR for low-precision training by increasing the batch size. This shrinks
the variance of the gradient distribution in Figure 1 without changing the mean and concentrates
more of the gradient distribution towards downhill directions, making the algorithm more greedy.
To verify this, we tried different batch sizes for SR including 128, 256, 512 and 1024, and found that
the larger the batch size, the better the performance of SR. Figure 5(a) illustrates the effect of a batch
size of 1024 for BC and SR methods. We find that the BC method, like classical SGD, performs best
8

30
20
10
0

0

20

40

60

80

100

120

Epochs

140

160

180

50

conv_1
conv_2
conv_3
conv_4
conv_5
conv_6
linear_1
linear_2
linear_3

40
30
20
10
0

0

20

(a) R-ADAM

40

60

80

100

120

Epochs

140

160

180

Percentage of changed weights (%)

conv_1
conv_2
conv_3
conv_4
conv_5
conv_6
linear_1
linear_2
linear_3

40

Percentage of changed weights (%)

Percentage of changed weights (%)

50

50

conv_1
conv_2
conv_3
conv_4
conv_5
conv_6
linear_1
linear_2
linear_3

40
30
20
10
0

0

20

(b) BC-ADAM

40

60

80

100

120

Epochs

140

160

180

(c) SR-ADAM

BC-ADAM 128
BC-ADAM 1024
SR-ADAM 128
SR-ADAM 1024

Error (%)

50
40
30
20
10
0

0

20

40

60

80

100

Epochs

120

140

160

(a) BC-ADAM vs SR-ADAM

60
50
40
30

BC-ADAM 128
BC-ADAM 1024
SR-ADAM 128
SR-ADAM 1024

20
10
0

0

20

40

60

80

100

Epochs

120

140

160

Percentage of changed weights (%)

60

Percentage of changed weights (%)

Figure 4: Percentage of weight changes during training of VGG-BC on CIFAR-10.
50

BC-ADAM 128
BC-ADAM 1024
SR-ADAM 128
SR-ADAM 1024

40
30
20
10
0

0

20

40

60

80

100

Epochs

120

140

160

(b) Weight changes since beginning (c) Weight changes every 5 epochs

Figure 5: Effect of batch size on SR-ADAM when tested with ResNet-56 on CIFAR-10. (a) Test error vs epoch.
Test error is reported with dashed lines, train error with solid lines. (b) Percentage of weight changes since
initialization. (c) Percentage of weight changes per every 5 epochs.
with a small batch size. However, a large batch size is essential for the SR method to perform well.
Figure 5(b) shows the percentage of weights changed by SR and BC during training. We see that the
large batch methods change the weights less aggressively than the small batch methods, indicating
less exploration. Figure 5(c) shows the percentage of weights changed during each 5 epochs of
training. It is clear that small-batch SR changes weights much more frequently than using a big batch.
This property of big batch training clearly benefits SR; we see in Figure 5(a) and Table 1 that big
batch training improved performance over SR-ADAM consistently.
In addition to providing a means of improving fixed-point training, this suggests that recently
proposed methods using big batches [28, 29] may be able to exploit lower levels of precision to
further accelerate training.

7

Conclusion

The training of quantized neural networks is essential for deploying machine learning models
on portable and ubiquitous devices. We provide a theoretical analysis to better understand the
BinaryConnect (BC) and Stochastic Rounding (SR) methods for training quantized networks. We
proved convergence results for BC and SR methods that predict an accuracy bound that depends
on the coarseness of discretization. For general non-convex problems, we proved that SR differs
from conventional stochastic methods in that it is unable to exploit greedy local search. Experiments
confirm these findings, and show that the mathematical properties of SR are indeed observable (and
very important) in practice.

Acknowledgments
T. Goldstein was supported in part by the US National Science Foundation (NSF) under grant CCF1535902, by the US Office of Naval Research under grant N00014-17-1-2078, and by the Sloan
Foundation. C. Studer was supported in part by Xilinx, Inc. and by the US NSF under grants
ECCS-1408006, CCF-1535897, and CAREER CCF-1652065. H. Samet was supported in part by the
US NSF under grant IIS-13-20791.
9

References
[1] Courbariaux, M., Hubara, I., Soudry, D., El-Yaniv, R., Bengio, Y.: Binarized neural networks: Training
deep neural networks with weights and activations constrained to +1 or -1. arXiv preprint arXiv:1602.02830
(2016)
[2] Marchesi, M., Orlandi, G., Piazza, F., Uncini, A.: Fast neural networks without multipliers. IEEE
Transactions on Neural Networks 4(1) (1993) 53–62
[3] Rastegari, M., Ordonez, V., Redmon, J., Farhadi, A.: XNOR-Net: ImageNet Classification Using Binary
Convolutional Neural Networks. ECCV (2016)
[4] Gupta, S., Agrawal, A., Gopalakrishnan, K., Narayanan, P.: Deep learning with limited numerical precision.
In: ICML. (2015)
[5] Courbariaux, M., Bengio, Y., David, J.P.: Binaryconnect: Training deep neural networks with binary
weights during propagations. In: NIPS. (2015)
[6] Lin, D., Talathi, S., Annapureddy, S.: Fixed point quantization of deep convolutional networks. In: ICML.
(2016)
[7] Hwang, K., Sung, W.: Fixed-point feedforward deep neural network design using weights+ 1, 0, and- 1.
In: IEEE Workshop on Signal Processing Systems (SiPS). (2014)
[8] Lin, Z., Courbariaux, M., Memisevic, R., Bengio, Y.: Neural networks with few multiplications. ICLR
(2016)
[9] Li, F., Zhang, B., Liu, B.: Ternary weight networks. arXiv preprint arXiv:1605.04711 (2016)
[10] Kim, M., Smaragdis, P.: Bitwise neural networks. In: ICML Workshop on Resource-Efficient Machine
Learning. (2015)
[11] Hubara, I., Courbariaux, M., Soudry, D., El-Yaniv, R., Bengio, Y.: Quantized neural networks: Training
neural networks with low precision weights and activations. arXiv preprint arXiv:1609.07061 (2016)
[12] Baldassi, C., Ingrosso, A., Lucibello, C., Saglietti, L., Zecchina, R.: Subdominant dense clusters allow for
simple learning and high computational performance in neural networks with discrete synapses. Physical
review letters 115(12) (2015) 128101
[13] Anwar, S., Hwang, K., Sung, W.: Fixed point optimization of deep convolutional neural networks for
object recognition. In: ICASSP, IEEE (2015)
[14] Zhu, C., Han, S., Mao, H., Dally, W.J.: Trained ternary quantization. ICLR (2017)
[15] Zhou, A., Yao, A., Guo, Y., Xu, L., Chen, Y.: Incremental network quantization: Towards lossless CNNs
with low-precision weights. ICLR (2017)
[16] Zhou, S., Wu, Y., Ni, Z., Zhou, X., Wen, H., Zou, Y.: Dorefa-net: Training low bitwidth convolutional
neural networks with low bitwidth gradients. arXiv preprint arXiv:1606.06160 (2016)
[17] Höhfeld, M., Fahlman, S.E.: Probabilistic rounding in neural network learning with limited precision.
Neurocomputing 4(6) (1992) 291–299
[18] Miyashita, D., Lee, E.H., Murmann, B.: Convolutional neural networks using logarithmic data representation. arXiv preprint arXiv:1603.01025 (2016)
[19] Soudry, D., Hubara, I., Meir, R.: Expectation backpropagation: Parameter-free training of multilayer
neural networks with continuous or discrete weights. In: NIPS. (2014)
[20] Cheng, Z., Soudry, D., Mao, Z., Lan, Z.: Training binary multilayer neural networks for image classification
using expectation backpropagation. arXiv preprint arXiv:1503.03562 (2015)
[21] Martens, J., Grosse, R.: Optimizing neural networks with kronecker-factored approximate curvature. In:
International Conference on Machine Learning. (2015) 2408–2417
[22] Lan, G., Nemirovski, A., Shapiro, A.: Validation analysis of mirror descent stochastic approximation
method. Mathematical programming 134(2) (2012) 425–458
[23] Levin, D.A., Peres, Y., Wilmer, E.L.: Markov chains and mixing times. American Mathematical Soc.
(2009)

10

[24] Simonyan, K., Zisserman, A.: Very Deep Convolutional Networks for Large-Scale Image Recognition. In:
ICLR. (2015)
[25] He, K., Zhang, X., Ren, S., Sun, J.: Deep Residual Learning for Image Recognition. In: CVPR. (2016)
[26] Russakovsky, O., Deng, J., Su, H., Krause, J., Satheesh, S., Ma, S., Huang, Z., Karpathy, A., Khosla, A.,
Bernstein, M., et al.: Imagenet Large Scale Visual Recognition Challenge. IJCV (2015)
[27] Kingma, D., Ba, J.: Adam: A method for stochastic optimization. ICLR (2015)
[28] De, S., Yadav, A., Jacobs, D., Goldstein, T.: Big batch SGD: Automated inference using adaptive batch
sizes. arXiv preprint arXiv:1610.05792 (2016)
[29] Goyal, P., Dollár, P., Girshick, R., Noordhuis, P., Wesolowski, L., Kyrola, A., Tulloch, A., Jia, Y., He, K.:
Accurate, large minibatch sgd: Training imagenet in 1 hour. arXiv preprint arXiv:1706.02677 (2017)
[30] Lax, P.: Linear Algebra and Its Applications. Number v. 10 in Linear algebra and its applications. Wiley
(2007)
[31] Krizhevsky, A.: Learning multiple layers of features from tiny images. (2009)
[32] Zagoruyko, S., Komodakis, N.: Wide residual networks. arXiv preprint arXiv:1605.07146 (2016)
[33] Collobert, R., Kavukcuoglu, K., Farabet, C.: Torch7: A matlab-like environment for machine learning. In:
BigLearn, NIPS Workshop. (2011)
[34] Ioffe, S., Szegedy, C.: Batch Normalization: Accelerating Deep Network Training by Reducing Internal
Covariate Shift. (2015)

11

Training Quantized Nets: A Deeper Understanding
Appendices
Here we present proofs of the lemmas and theorems presented in the main paper, as well as some additional
experimental details and results.

A

Proof of Lemma 1

Proof. We want to bound the quantization error rt . Consider the i-th entry in rt denoted by rit . Similarly, we
define wit and ∇f˜(wt )i . Choose some random number p ∈ [0, 1]. The stochastic rounding operation produces a
value of rt given by
rit = Qs (wit − αt ∇f˜(wt )i ) − wit + αt ∇f˜(wt )i





αt ∇f˜(wt )i
−αt ∇f˜(wt )i
−αt ∇f˜(wt )i
αt ∇f˜(wt )i


+
−
+
1,
for
p
≤
−
,


∆
∆
∆
∆
=∆·
 αt ∇f˜(wt )i  −αt ∇f˜(wt )i 



+
,
otherwise,
∆
∆
(
−q + 1, for p ≤ q,
=∆·
−q,
otherwise,


−αt ∇f˜(wt )i
αt ∇f˜(wt )i
−
and q ∈ [0, 1].
where we write q = −
∆
∆
Now we have


Ep (rit )2 ≤ ∆2 ((−q + 1)2 q + (−q)2 (1 − q))
= ∆2 q(1 − q)
≤ ∆2 min{q, 1 − q}.
Because min{q, 1 − q} ≤



αt ∇f˜(wt )i
αt ∇f˜(wt )i
, it follows that Ep (rit )2 ≤ ∆2
≤ ∆ αt ∇f˜(wt )i .
∆
∆

Summing over the index i yields
Ep rt

Now, E ∇f˜(wt )

B

2
2

≤ E ∇f˜(wt )

2
2

2
2

≤ ∆αt ∇f˜(wt ) 1
√
≤ dαt ∆ ∇f˜(wt )

2

.

≤ G2 . Plugging this into (11) yields
√
2
E rt 2 ≤ d∆αt G.

(11)

(12)

Proof of Theorem 1

Proof. From the update rule (6), we get:

wt+1 = Q wt − αt ∇f˜(wt )
= wt − αt ∇f˜(wt ) + rt ,
where rt denotes the quantization used on the t-th iteration. Subtracting by the optimal w? , taking norm, and
taking expectation conditioned on wt , we get:
Ekwt+1 − w? k2 = kwt − w? k2 − 2Ehwt − w? , αt ∇f˜(wt ) − rt i + Ekαt ∇f˜(wt ) − rt k2
= kwt − w? k2 − 2αt hwt − w? , ∇F (wt )i + αt2 Ek∇f˜(wt )k2 + Ekrt k2
√
≤ kwt − w? k2 − 2αt hwt − w? , ∇F (wt )i + αt2 G2 + d∆αt G,
where we use the bounded variance assumption, E[rt ] = 0, and Lemma 1. Using the assumption that F is
µ-strongly convex, we can simplify this to:
√
Ekwt+1 − w? k2 ≤ (1 − αt µ)kwt − w? k2 − 2αt (F (wt ) − F (w? )) + αt2 G2 + d∆αt G.

12

Re-arranging the terms, and taking expectation we get:

√
2αt E(F (wt ) − F (w? )) ≤ (1 − αt µ)Ekwt − w? k2 − Ekwt+1 − w? k2 + αt2 G2 + d∆αt G.
√


1
µ
αt 2
d∆G
1
−
Ekwt+1 − w? k2 +
G +
.
⇒ E(F (wt ) − F (w? )) ≤
Ekwt − w? k2 −
2αt
2
2αt
2
2

Assume that the stepsize decreases with the rate αt = 1/µ(t + 1). Then we have:
E(F (wt ) − F (w? )) ≤

µ(t + 1)
µt
1
Ekwt − w? k2 −
Ekwt+1 − w? k2 +
G2 +
2
2
2µ(t + 1)

√

d∆G
.
2

Averaging over t = 0 to T , we get a telescoping sum on the right hand side, which yields:
√
T
T
µ(T + 1)
1 X
G2 X 1
d∆G
t
?
+
−
EkwT +1 − w? k2
E(F (w ) − F (w )) ≤
T t=0
2µT t=0 t + 1
2
2
√
(1 + log(T + 1))G2
d∆G
≤
+
.
2µT
2
Using Jensen’s inequality, we have:
E(F (w̄T ) − F (w? )) ≤
where w̄T =

1
T

PT

t=0

T
1 X
E(F (wt ) − F (w? )),
T t=0

wt , the average of the iterates.

Thus the final convergence theorem is given by:
(1 + log(T + 1))G2
+
E[F (w̄ ) − F (w )] ≤
2µT
T

C

√

?

d∆G
.
2

Proof of Theorem 2

Proof. From the update rule (6), we have,

wt+1 = Q wt − αt ∇f˜(wt ) = wt − αt ∇f˜(wt ) + rt ,
where rt denotes the quantization error on the t-th iteration. Hence we have
kwt+1 − w? k2 = kwt − αt ∇f˜(wt ) + rt − w? k2
= kwt − w? k2 − 2hwt − w? , αt ∇f˜(wt ) − rt i + kαt ∇f˜(wt ) − rt k2 .
Taking expectation, and using E[f˜(wt )] = ∇F (wt ) and E[rt ] = 0, we have
Ekwt+1 − w? k2 = Ekwt − w? k2 − 2αt Ehwt − w? , ∇F (wt )i + Ekαt ∇f˜(wt ) − rt k2
= Ekwt − w? k2 − 2αt Ehwt − w? , ∇F (wt )i + αt2 Ek∇f˜(wt )k2 + Ekrt k2 .
Using the bounded variance assumption Ek∇f˜(wt )k2 ≤ G2 and bounded quantization error in Lemma 1, we
have
√
Ekwt+1 − w? k2 ≤ Ekwt − w? k2 − 2αt Ehwt − w? , ∇F (wt )i + αt2 G2 + d∆αt G.
(13)
F (x) is convex and hence h∇F (x), xt − x? i ≥ F (xt ) − F (x∗ ), which can be used in (13) to get
√
Ekwt+1 − w? k2 ≤ Ekwt − w? k2 − 2αt E[F (wt ) − F (w? )] + αt2 G2 + d∆αt G.
Re-arranging the terms, we have,
 αt 2 1 √
1
Ekwt − w? k2 − Ekwt+1 − w? k2 +
G +
d∆G.
2αt
2
2
Accumulate from t = 1 to T to get

T 
T
X
X
1
1
1
Ekw1 − w? k2 +
−
E[F (wt ) − F (w? )] ≤
Ekwt − w? k2
2α1
2αt
2αt−1
t=1
t=1
E[F (wt ) − F (w? )] ≤

+

T
X
αt 2 T √
G +
d∆G.
2
2
t=1

13

√
αt ≤ c T + 1, we have
√
√
T
X
T 2 c T +1 2 T√
t
?
D +
G +
d∆G.
E[F (w ) − F (w )] ≤
2c
2
2
t=1

Applying Ekwt − w? k2 ≤ D2 and

PT

t=1

Since F (w) is convex, we can set w̄T =

1
T

PT

t=1

(14)

wt , and use Jensen’s inequality to arrive at

E[F (w̄T ) − F (w? )] ≤

T
1 X
E[F (wt ) − F (w? )].
T t=1

(15)

Combine (14) and (15) to achieve
1
E[F (w̄ ) − F (w )] ≤ √ D2 +
2c T
T

D

?

√
√
T +1 2
d∆G
cG +
.
2T
2

Proof of Theorem 3

Proof. From the update rule (7), we have,


wt+1 = wt − αt ∇f˜ Q(wt ) = wt − αt ∇f˜ wt + rt .
Taking expectation conditioned on wt and rt , we have
Ekwt+1 − w? k2

= Ekwt − αt ∇f˜ wt + rt − w? k2


= Ekwt − αt ∇f˜ wt ) + αt ∇f˜ wt − αt ∇f˜ wt + rt − w? k2




= kwt − w? k2 − 2αt Ehwt − w? , ∇f˜ wt i + 2αt Ehwt − w? , ∇f˜ wt − ∇f˜ wt + rt i + Ekαt ∇f˜ wt + rt k2




= kwt − w? k2 − 2αt hwt − w? , ∇F wt i + 2αt hwt − w? , ∇F wt − ∇F wt + rt i + αt2 Ek∇f˜ wt + rt k2



≤ kwt − w? k2 − 2αt hwt − w? , ∇F wt i + 2αt kwt − w? kk∇F wt − ∇F wt + rt k + αt2 G2

≤ kwt − w? k2 − 2αt hwt − w? , ∇F wt i + 2αt Lkrt kkwt − w? k + αt2 G2 .
√
Using krt k ≤ d∆ and the bounded domain assumption, we get
√

Ekwt+1 − w? k2 ≤ kwt − w? k2 − 2αt hwt − w? , ∇F wt i + 2αt L d∆kwt − w? k + αt2 G2
√

≤ kwt − w? k2 − 2αt hwt − w? , ∇F wt i + 2αt L d∆D + αt2 G2 .
Taking expectation, and following the same steps as in Theorem 2, we get the convergence result:
√
1
T +1 2 √
E[F (w̄T ) − F (w? )] ≤ √ D2 +
cG + d∆LD.
2T
2c T

E

Proof of Theorem 4

Proof. From the update rule (7), we get

wt+1 = wt − αt ∇f˜ Q(wt )

= wt − αt ∇f˜ wt + rt


= wt − αt [∇f˜ wt + ∇2 f˜ wt rt + r̂t ]
where kr̂t k ≤ L22 krt k2 from our assumption on the Hessian. Note that in general rt has mean zero while r̂t
does not. Using the same steps as in the Theorem 1, we get
Ekwt+1 − w? k2 = kwt − w? k2 − 2αt Ehwt − w? , ∇f˜(wt + rt )i + αt2 Ek∇f˜(wt + rt )k2 .
≤ kwt − w? k2 − 2αt Ehwt − w? , ∇F (wt ) + r̂t i + αt2 G2
= kwt − w? k2 − 2αt Ehwt − w? , ∇F (wt )i + αt2 G2 − 2αt Ehwt − w? , r̂t i

14

Assuming the domain has finite
√ diameter D, and observing that the quantization error for BC-SGD can always
be upper-bounded as krt k ≤ d∆, we get:
√
L2
−2αt Ehwt − w? , r̂t i ≤ 2αt DEkr̂t k ≤ 2αt D krt k ≤ αt DL2 d∆.
2
Following the same steps as in Theorem 1, we get
√
(1 + log(T + 1))G2
DL2 d∆
T
?
+
.
E[F (w̄ ) − F (w )] ≤
2µT
2

F

Proof of Theorem 5

Proof. Let the matrix Uα be a partial transition matrix defined by Uα (x, x) = 0, and Uα (x, y) = Tα (x, y) for
x 6= y. From Uα , we can get back the full transition matrix Tα using the formula
Tα = I − diag(1T Uα ) + Uα .
Note that this formula is essentially “filling in” the diagonal entries of Tα so that every column sums to 1, thus
making Tα a valid stochastic matrix.
Let’s bound the entries in Uα . Suppose that we begin an iteration of the stochastic rounding algorithm at some
point x. Consider an adjacent point y that differs from x at only 1 coordinate, k, with yk = xk + ∆. Then we
have

Uα (x, y) =

1
α

∆

Z

px,k (x/α)
0

1
x
dx +
∆
α

Z

2∆

px,k (x/α)
∆

2∆ − x
dx
∆

Z
αz
1 2∆/α
2∆ − αz
px,k (z) α dz +
px,k (z)
α dz
∆
α ∆/α
∆
0
Z ∞
Z ∆/α
z
≤α
px,k (z) dz
px,k (z) dz +
∆
∆/α
0
Z ∞
z
=α
px,k (z) dz + O(α2 ).
∆
0
Note we have used the decay assumption:
Z ∞
C
px,k (z) ≤ 2 .
ν
ν
Likewise, if yk = xk − ∆, then the transition probability is
Z 0
z
Uα (x, y) = α
px,k (z) dz + O(α2 ),
∆
−∞
1
=
α

Z

∆/α

(16)

(17)

and if yk = xk ± m∆ for an integer m > 1,
Uα (x, y) = O(α2 ).

(18)

We can approximate the behavior of Uα using the matrix
R ∞
z

R0 px,k (z) ∆ dz, if x and y differ only at coordinate k, and yk = xk + ∆
0
z
Ũ (x, y) =
p (z) ∆
dz, if x and y differ only at coordinate k, and yk = xk − ∆
−∞ x,k

0, otherwise.
Define the associated markov chain transition matrix
T̃α0 = I − α0 · diag(1T Ũ ) + α0 Ũ ,

(19)

where α0 is the largest scalar such that the stochastic linear operator T̃α0 has non-negative entries. For α < α0 ,
T̃α has non-negative entries and column sums equal to 1; it thus defines the transition operator of a markov
chain. Let π̃ denote the stationary distribution of the markov chain with transition matrix T̃α0 .
We now claim that π̃ is also the stationary distribution of T̃α for all α < α0 . We verify this by noting that
T̃α = (I − α · diag(1T Ũ )) + αŨ
α
α
= (1 −
)I +
[I − α0 · diag(1T Ũ ) + α0 Ũ ]
α0
α0
α
α
= (1 −
)I +
T̃α ,
α0
α0 0

15

(20)

and so T̃α π̃ = (1 −

α
)π̃
α0

+

α
π̃
α0

= π̃.

Recall that Tα is the transition matrix for the Markov chain generated by the stochastic rounding algorithm with
learning rate α. We wish to show that this markov chain is well approximated by T̃α . Note that
Y
Tα (x, y) =
Tα (x, x + (yk − xk )∆ek ) ≤ O(α2 )
k,xk 6=yk

when x, y differ at more than 1 coordinate. In other words, transitions between multiple coordinates simultaneously become vanishingly unlikely for small α. When x and y differ by exactly 1 coordinate, we know from
(16) that
Tα (x, y) = αU (x, y) + O(α2 ).
These observations show that the off-diagonal elements of Tα are well approximated (up to uniform O(α2 )
error) by the corresponding elements in αU. Since the columns of Tα sum to one, the diagonal elements are well
approximated as well, and we have
Tα = (I − α · diag(1T U )) + αU + O(α2 ) = T̃α + O(α2 ).
To be precise, the notation above means that
|Tα (x, y) − T̃α (x, y)| < Cα2 ,

(21)

for some C that is uniform over (x, y).
We are now ready to show that the stationary distribution of Tα exists and approaches π̃. Re-arranging (20) gives
us
α0 T̃α + (α − α0 )I = αT̃α0 .
Combining this with (21), we get
α0 Tα + (α − α0 )I − αT̃α0

∞

α0
α0
Tα + (1 −
)I − T̃α0
α
α

∞

< O(α2 ),

and so
< O(α).

(22)

From (22), we see that the matrix αα0 Tα + (1 − αα0 )I approaches T̃α0 . Note that π̃ is the Perron-Frobenius
eigenvalue of T̃α0 , and thus has multiplicity 1. Multiplicity 1 eigenvalues/vectors of a matrix vary continuously
with small perturbations to that matrix (Theorem 8, p130 of [30]). It follows that, for small α, αα0 Tα + (1 − αα0 )I
has a stationary distribution, and this distribution approaches π̃. The leading eigenvector of αα0 Tα + (1 − αα0 )I is
the same as the leading eigenvector of Tα , and it follows that Tα has a stationary distribution that approaches π̃.
RC
R0
Finally, note that we have assumed 0 2 px,k (z) dz > 0 and −C2 px,k (z) dz > 0. Under this assumption, for
α<

1
,
C2

T̃α0 (x, y) > 0 whenever x, y are neighbors the differ at a single coordinate. It follows that every

state in the Markov chain T̃α0 is accessible from every other state by traversing a path of non-zero transition
probabilities, and so π̃(x) > 0 for every state x.

G

Proof of Theorem 6

Proof. Given some distribution π over the states of the markov chain, and some set A of states, let [π]A =
P
a∈A π(a) denote the measure of A with respect to π.
Suppose for contradiction that the mixing time of the chain remains bounded as α vanishes. Then we can find an
integer M that upper bounds the -mixing time for all α. By the assumption of the theorem, we can select some
set of states A with [π̃]A > , and some starting state y 6∈ A. Let e be a distribution (a vector in the finite-state
case) with ey = 1, ek = 0 for k 6= y. Note that [e]A = 0 because y 6∈ A. Then
[e]A − [π̃]A > .
Note that, as α → 0, we have kTα − T̃α k → 0 and thus kTαM − T̃αM k → 0. We also see from the definition
of T̃α in (19), limα→0 T̃α = I. It follows that
lim [TαM e]A − [π̃]A = [e]A − [π̃]A > ,

α→0

and so for some α the inequality (10) is violated. This is a contradiction because it was assumed M is an upper
bound on the mixing time.

16

H
H.1

Additional Experimental Details & Results
Neural Net Architecture & Training Details

We train image classifiers using two types of networks, VGG-like CNNs [24] and Residual networks [25], on
CIFAR-10/100 [31] and ImageNet 2012 [26]. VGG-9 on CIFAR-10 consists of 7 convolutional layers and 2 fully
connected layers. The convolutional layers contain 64, 64, 128, 128, 256, 256 and 256 of 3 × 3 filters respectively.
There is a Batch Normalization and ReLU after each convolutional layer and the first fully connected layer. The
details of the architecture are presented in Table 2. VGG-BC is a high-capacity network used for the original BC
method [5], which contains 6 convolutional layers and 3 linear layers. We use the same architecture as in [5]
except using softmax and cross-entropy loss instead of SVM and squared hinge loss, respectively. The details
of the architecture are presented in Table 3. ResNets-56 has 55 convolutional layers and one linear layer, and
contains three stages of residual blocks where each stage has the same number of residual blocks. We also create
a wide ResNet-56 (WRN-56-2) that doubles the number of filters in each residual block as in [32]. ResNets-18
for ImageNet has the same description as in [25].
The default minibatch size is 128. However, the big-batch SR-ADAM method adopts a large minibatch size (512
for WRN-56-2 and ResNet-18 and 1024 for other models). Following [5], we do not use weight decay during
training. We implement all models in Torch7 [33] and train the quantized models with NVIDIA GPUs.
Similar to [3], we only quantize the weights in the convolutional layers, but not linear layers, during training.
Binarizing linear layers causes some performance drop without much computational speedup. This is because
fully connected layers have very little computation overhead compared to Conv layers. Also, for state-of-the-art
CNNs, the number of FC parameters is quite small. The number of params of Conv/FC layers for CNNs in Table
1 are (in millions): VGG-9: 1.7/1.1, VGG-BC: 4.6/9.4, ResNet-56: 0.84/0.0006, WRN-56-2: 3.4/0.001, ResNet18: 11.2/0.5. While the VGG-like nets have many FC parameters, the more efficient and higher performing
ResNets are almost entirely convolutional.

Table 2: VGG-9 on CIFAR-10.
layer type
kernel size
Conv_1
3×3
Conv_2
3×3
Max Pooling
2×2
Conv_3
3×3
Conv_4
3×3
2×2
Max Pooling
Conv_5
3×3
Conv_6
3×3
3×3
Conv_7
Max Pooling
2×2
Linear
1×1
1×1
Linear

input size
3 × 32 × 32
64 × 32 × 32
64 × 32 × 32
64 × 16 × 16
128 × 16 × 16
128 × 16 × 16
128 × 8 × 8
256 × 8 × 8
256 × 8 × 8
256 × 8 × 8
1 × 4096
1 × 256

output size
64 × 32 × 32
64 × 32 × 32
64 × 16 × 16
128 × 16 × 16
128 × 16 × 16
128 × 8 × 8
256 × 8 × 8
256 × 8 × 8
256 × 8 × 8
256 × 4 × 4
1 × 256
1 × 10

Table 3: VGG-BC for CIFAR-10.
layer type
kernel size
Conv_1
3×3
3×3
Conv_2
Max Pooling
2×2
Conv_3
3×3
Conv_4
3×3
Max Pooling
2×2
Conv_5
3×3
Conv_6
3×3
Max Pooling
2×2
Linear
1×1
Linear
1×1
Linear
1×1

H.2

input size
3 × 32 × 32
128 × 32 × 32
128 × 32 × 32
128 × 16 × 16
256 × 16 × 16
256 × 16 × 16
256 × 8 × 8
512 × 8 × 8
512 × 8 × 8
1 × 8192
1 × 1024
1 × 1024

output size
128 × 32 × 32
128 × 32 × 32
128 × 16 × 16
256 × 16 × 16
256 × 16 × 16
256 × 8 × 8
512 × 8 × 8
512 × 8 × 8
512 × 4 × 4
1 × 1024
1 × 1024
1 × 10

Convergence Curves

The convergence curves for training and testing errors reported in Table 1 are shown in Figure 6.

17

30
20

40

10

0

20

40

60

80

120

140

20

40

0

160

ADAM 128
R-ADAM 128
BC-ADAM 128
SR-ADAM 512

50

Error (%)

100

Epochs
(a) VGG-9 on CIFAR-10

60

30
20

0

20

40

60

80

100

120

140

Epochs
(b) VGG-BC on CIFAR-10

40

60

80

100

120

140

80
70
60

Epochs
(d) WSN-56-2 on CIFAR-10

0

20

40

60

80

100

140

160

ADAM 128
R-ADAM 128
BC-ADAM 128
SR-ADAM 512

80

40

120

Epochs
(c) ResNet-56 on CIFAR-10

90

50

10

160

20

100

ADAM 128
R-ADAM 128
BC-ADAM 128
SR-ADAM 1024

90

70
60
50
40

20

20

30

0

160

30

0

40

10

100

10
0

30

ADAM 128
R-ADAM 128
BC-ADAM 128
SR-ADAM 1024

50

10

Error (%)

0

60

ADAM 128
R-ADAM 128
BC-ADAM 128
SR-ADAM 1024

Error (%)

Error (%)

40

50

Error (%)

ADAM 128
R-ADAM 128
BC-ADAM 128
SR-ADAM 1024

Error (%)

50

0

20

40

60

80

100

120

140

Epochs
(e) ResNet-56 on CIFAR-100

160

30

0

10

20

30

40

50

60

70

80

90

Epochs
(f) ResNet-18 on ImageNet 2012

Figure 6: Training and testing errors of different training methods for VGG-9, VGG-BC, ResNet-56, WideResNet-56-2 and ResNet-18. The solid line is the training error and the dashed line is the testing error.

H.3

Weight Initialization and Learning Rate

For experiments on SR-Adam and R-Adam, the weights of convolutional layers are intitialized with random
Rademacher (±1) variables. The authors of BC [5] adopt a small initial learning rate (0.003) and it takes 500
epochs to converge. It is observed that large binary weights (∆ = 1) will generate small gradients when batch
normalization is used [34], hence a large learning rate is necessary for faster convergence. We experiment with a
larger learning rate (0.01) and find it converges to the same performance within 160 epochs, comparing with 500
epochs in the original paper [5].

H.4

Weight Decay

Figure 7 shows the effect of applying weight decay to BC-ADAM. As shown in Figure 7(a), BC-ADAM
with 1e-5 weight decay yields worse performance compared to zero weight decay. Applying weight decay in
BC-ADAM will shrink wr to 0, as well as increase the distance between wb and wr . Figure 7(b) and 7(c) shows
the distance between wb and wr during training. With 1e-5 weight decay, the average weight difference between
wb and wr approaches 1, which indicates wr is close to zero. Weight decay cannot “decay” the weight of SR as
kwb k2 is the same for all binarized networks.

Error (%)

40
30
20
10
0

0

20

40

60

80

100

120

140

Epochs
(a) WD=1e-5 vs WD=0

160

1.4

2.5

1.2

conv_1
conv_2
conv_3
conv_4
conv_5
conv_6
linear_1
linear_2
linear_3

1.0
0.8
0.6
0.4
0.2
0.0

0

20

(b)

40

|wbt

60

−

80

100

Epochs

wrt |,

120

140

160

WD=1e-5

Average weight distance

BC-ADAM 128 WD=1e-5
BC-ADAM 128 WD=0

Average weight distance

50

180

conv_1
conv_2
conv_3
conv_4
conv_5
conv_6
linear_1
linear_2
linear_3

2.0
1.5
1.0
0.5
0.0

0

20

(c)

40

60

|wbt

80

100

Epochs

−

wrt |,

120

140

160

WD=0

Figure 7: The effect of weight decay (WD) on BC-ADAM for training VGG-BC. The y-axis of (b) and (c) is the

averaged weight difference between the binary weights wb and the real-valued weights wr , i.e., d1 kwbt − wrt k1 .
where d is the number of weights in wb .

18

180

