Batch Normalization: Accelerating Deep Network Training by Reducing Internal Covariate Shift

Standardize every layer's inputs over each mini-batch, and an ImageNet classifier reaches its original accuracy in 14 times fewer training steps.

Ioffe and Szegedy's 2015 layer subtracts the batch mean, divides by the batch standard deviation, then applies a learned scale and shift. It lets networks train at learning rates 5 to 30 times higher, with less careful weight initialization, and makes sigmoid networks trainable again, though the paper's explanation for why it works has not held up.

Four lines of arithmetic per feature, the chain rule through a batch average, a second set of statistics for test time, and an argument about weight scale for why the layer tolerates high learning rates.

A deep network is trained by stochastic gradient descent (SGD): take a mini-batch of examples, compute the loss, compute its gradient with respect to every weight, move each weight a small step against its gradient, and repeat for millions of steps. Before 2015, deep networks trained this way needed small learning rates and careful weight initialization, and networks with saturating nonlinearities such as the sigmoid were close to untrainable at depth.

Ioffe and Szegedy attributed much of this to one effect. A layer's input is the previous layer's output, so every update to an early layer changes the distribution of values arriving at every later layer. They named this internal covariate shift and proposed a layer that holds the mean and variance of each input fixed during training. Santurkar et al. (2018) later tested that explanation and found it does not account for the speedup; the last section covers their experiments.

The page runs in dependency order: why a drifting input distribution stalls a sigmoid layer, why normalizing outside the gradient computation fails, the four-line transform itself, how gradients flow through the batch mean and variance, what changes at inference time, why the layer tolerates much higher learning rates, the ImageNet numbers, and what later papers found about the mechanism.

Why deep networks trained slowly: saturation and shifting inputs

The sigmoid $g(x) = 1/(1+e^{-x})$ maps any real number into $(0,1)$. Its derivative is:

which is $0.25$ at $x = 0$, $0.0177$ at $x = 4$ and $0.0066$ at $x = 5$. Backpropagation multiplies the gradient by this slope on the way down, so an example whose pre-activation sits at $x = 5$ passes back 38 times less gradient than one at $x = 0$. If most of a batch sits out there, the weights feeding that unit barely change. This is saturation, and the paper's Section 2 lists the usual remedies of the time: ReLU units, careful initialization, and small learning rates.

The pre-activation is $x = Wu + b$, where $u$ is the output of the layer below. When the layers below change, the distribution of $u$ changes, so a unit that started in the steep middle of the sigmoid can drift into a flat tail during training. The paper calls a change in a model's input distribution covariate shift, a term from Shimodaira (2000), and defines internal covariate shift as "the change in the distribution of network activations due to the change in network parameters during training".

Section 4.1 of the paper shows the effect on MNIST with a small network: 3 fully connected hidden layers of 100 sigmoid units each, weights initialized to small random Gaussian values, 50,000 steps at 60 examples per mini-batch. For one typical unit in the last hidden layer, the 15th, 50th and 85th percentiles of its sigmoid input changed significantly in both mean and variance over training without BN and were much more stable with it, and the BN network reached higher test accuracy.

In the figure below, the drift slider moves the mean of a batch of pre-activations. Watch the gradient-signal readout as the batch moves into a tail, then switch Batch Norm on.

With Batch Norm on, the readout stays near its maximum at every drift, because the normalized batch always has mean 0 and variance 1 before it reaches the sigmoid. The paper's sigmoid ImageNet network depends on this: with BN it reaches 69.8% accuracy, and without BN it stays at chance (the ImageNet section has the details).

Why you can't just normalize outside the gradient

An obvious implementation computes each layer's mean from the data, separately from training, and subtracts it before each step, treating the mean as a constant. Section 2 of the paper shows why this fails, and the authors report that in early experiments the model blew up when the normalization statistics were computed outside the gradient step.

Take a layer that adds a learned bias $b$, so $x = u + b$, and normalizes by subtracting the mean over the training data, $\hat{x} = x - \mathrm{E}[x]$. If the gradient step treats $\mathrm{E}[x]$ as a constant, the update to $b$ is $\Delta b$ with $\Delta b \propto -\partial\ell/\partial\hat{x}$, the gradient of the loss with respect to the normalized value. After the update, the normalized output is:

Adding $\Delta b$ to every $x$ raises their mean by exactly $\Delta b$, and the subtraction removes it. The output is unchanged, so the loss is unchanged, so $\partial\ell/\partial\hat{x}$ is the same at the next step, and the next update moves $b$ by the same $\Delta b$ again. $b$ grows linearly, forever, while the loss stays fixed.

Scrub the training steps in the figure. In naive mode $b$ climbs while the layer output and the loss stay flat. In aware mode the gradient includes the dependence of the mean on $b$, which is exactly equal and opposite, so the gradient on $b$ is 0 and it stays put.

The paper concludes that normalization has to be part of the network, so that backpropagation differentiates through the mean and variance as functions of the parameters. In the paper's notation, the gradient must include $\partial\,\mathrm{Norm}/\partial X$ (the dependence on the whole batch) as well as $\partial\,\mathrm{Norm}/\partial x$. It also means the bias $b$ becomes useless inside a normalized layer, since any constant added before the mean subtraction is removed by it.

How batch normalization works

Batch Normalization treats each feature (one unit's pre-activation) separately. For one feature, a mini-batch $\mathcal{B}$ of $m$ examples gives $m$ values $x_1, \dots, x_m$. Algorithm 1 of the paper computes:

The variance divides by $m$, not $m-1$; this biased estimate is the one used in training, and the inference section returns to it. The constant $\epsilon$ is added to the variance inside the square root so the division is safe when a feature is nearly constant across the batch (PyTorch's default is $10^{-5}$). $\gamma$ and $\beta$ are learned, one pair per feature.

A worked example with $m = 4$, batch values 1, 2, 4 and 5, $\gamma = 2$, $\beta = 0.5$ and $\epsilon = 10^{-5}$:

The normalized values have mean 0 and variance 1 (up to $\epsilon$); the outputs have mean 0.5 and standard deviation 2.

Why not stop at $\hat{x}$? Forcing every sigmoid input to mean 0 and variance 1 would keep it in the sigmoid's near-linear middle and remove the network's ability to use the saturating parts when they help. The paper requires that the inserted transform can represent the identity, and $\gamma, \beta$ make it so:

exactly if $\epsilon$ is neglected. In the worked example, $\gamma = \sqrt{2.5} \approx 1.581$ and $\beta = 3$ give back $(1, 2, 4, 5)$. The network can therefore represent everything it could before, and gradient descent decides how far from the identity to go. Drag $\gamma$ and $\beta$ in the figure, or press "set to identity": the middle row never changes, and the bottom row can be placed anywhere.

Equation (2) is a simplified version of whitening, the classical input preprocessing that LeCun et al. (1998) and Wiesler and Ney (2011) found speeds up convergence: transform the inputs to zero mean and unit variance and also decorrelate them. Full whitening of a layer needs the covariance matrix of all its features and that matrix's inverse square root at every step, which is expensive and not everywhere differentiable; and with fewer examples in the batch than features, the batch covariance matrix is singular. The paper makes two simplifications. Each feature is normalized on its own, with no decorrelation. And the statistics come from the current mini-batch rather than the whole training set, so they are functions of the parameters that backpropagation can differentiate, which Figure 2 showed is required.

Where the layer goes, and convolutions

The paper inserts BN immediately before the nonlinearity. A layer $z = g(Wu + b)$ becomes:

with the bias $b$ dropped, because the mean subtraction cancels it and $\beta$ takes over its role. The paper normalizes $Wu$ rather than the layer input $u$: $u$ is the output of another nonlinearity, so its shape keeps changing, while $Wu$ is a sum of many terms and is "more likely to have a symmetric, non-sparse distribution, that is 'more Gaussian'", for which fixing the mean and variance fixes most of the distribution.

In a convolutional layer, a feature map applies the same filter at every position, and BN keeps that property by normalizing every position of a feature map the same way. The statistics for one feature map are pooled over the batch and all $p \times q$ positions, an effective batch of $m' = m\cdot p\cdot q$ values, and there is one $(\gamma, \beta)$ pair per feature map. At the paper's batch size of 32 and a $14 \times 14$ feature map (the size at the fourth stage of its Inception network), each mean and variance is computed over $32 \cdot 14 \cdot 14 = 6{,}272$ values. A convolution with 64 output channels gains 128 parameters, and those two per channel are all BN adds.

Backpropagating through the batch mean and variance

Every quantity in (2) is differentiable, including $\mu_{\mathcal B}$ and $\sigma^2_{\mathcal B}$. Because they are shared by the whole batch, example $i$'s output depends on every other example's input, and the gradient with respect to $x_i$ has three paths: the direct one through $\hat{x}_i$, one through the mean, and one through the variance. The paper gives the chain rule "before simplification":

To simplify (4), write $s$ for $\sqrt{\sigma^2_{\mathcal B} + \epsilon}$ and $g_i$ for $\partial\ell/\partial\hat{x}_i$. The second term of $\partial\ell/\partial\mu_{\mathcal B}$ contains the sum of the deviations $x_i - \mu_{\mathcal B}$, which is 0 by the definition of the mean, so only its first term survives. In the variance line, each deviation $x_i - \mu_{\mathcal B}$ equals $s\,\hat{x}_i$. Putting these into the last line of (4), the three paths contribute:

and together:

which is exact, $\epsilon$ included. Summing (5) over the batch, and summing it weighted by $\hat{x}_i$, gives two properties (the first uses that the $\hat{x}_i$ sum to 0, the second that their squares sum to $m\,\sigma^2_{\mathcal B}/(\sigma^2_{\mathcal B}+\epsilon)$):

Adding the same constant to every $x_i$, which is what a bias does, therefore has zero gradient (the term the naive scheme of Figure 2 left out), and stretching the batch about its mean has almost none.

Continuing the worked example (batch 1, 2, 4, 5 and $\gamma = 2$), suppose the loss gradient on the outputs is 1 for the first example and 0 for the others:

The four values sum to 0, and a finite-difference check on the forward pass gives the same numbers. Here the loss is $y_1$ itself. Raising $x_2$ or $x_3$ lowers $y_1$, mostly by raising the batch mean. Raising $x_4$ raises the mean too (a contribution of $-0.316$) but also widens the batch, which pulls the negative $\hat{x}_1$ toward 0 (a contribution of $+0.506$), so its net gradient is positive. If the loss gradient is 1 on every output, every $\partial\ell/\partial x_i$ is 0.

In the figure, the incoming gradient is a fixed random pattern plus a shared shift $c$ and a scale-shaped part $k\,\hat{x}_i$. With all three paths on, move either slider: the top bars change, the bottom bars do not. Cycle the paths button to direct-only and the shift passes straight through.

The gradients for $\gamma$ and $\beta$ are the plain sums in the last line of (4). Since the whole layer is a handful of means and elementwise operations, its backward pass costs about as much as its forward pass, and BN works with any gradient method; the paper names SGD with $m > 1$, momentum and Adagrad.

Training mode and inference mode

During training, $\mathrm{BN}(x_i)$ depends on which other examples share the batch. At inference the paper wants "the output to depend only on the input, deterministically", and a deployed model may see one example at a time. So after training, Algorithm 2 replaces the batch statistics with population statistics:

The expectations are averages over many training mini-batches of size $m$. The factor $m/(m-1)$ is Bessel's correction: a ÷$m$ batch variance measures spread around the batch's own mean, which is closer to the batch's points than the true mean is, so on average it comes out low by a factor $(m-1)/m$. For the worked batch, $2.5 \times 4/3 = 3.33$. At $m = 32$ the correction is $32/31$, about 3%; at $m = 2$ it doubles the estimate. The asymmetry is deliberate: training uses the ÷$m$ statistic that (2) defines and backprop differentiates, and inference uses the best estimate of the population variance.

The figure draws 40 random batches of size $m$ from a feature with mean 3 and variance 4. The top strip shows each batch's ÷$m$ variance and their average; the bottom strip follows one example, $x = 5$, through 40 different batches and shows the normalized value it gets in each. Slide $m$ down to 1 and up to 64.

Three things are visible. At $m = 4$ the average variance sits near 3, as $(m-1)/m \cdot 4$ predicts, and the corrected line near 4. At $m = 2$ the example's normalized value is always $\pm 1$, because any two distinct numbers standardize to $-1$ and $+1$; at $m = 1$ the batch variance is 0, every $\hat{x}$ is 0, the input is erased, and $m/(m-1)$ is undefined, which is why the paper requires $m > 1$. And at every batch size the training values scatter around the inference value 1.00, more widely for small $m$.

That scatter has two consequences. Section 3.4 of the paper treats it as a regularizer: an example is normalized differently each time it appears with different companions, and the authors found Dropout "can be either removed or reduced in strength". Their ImageNet setup also shuffled training data more thoroughly within shards so the same examples stopped landing in a batch together, worth about 1% validation accuracy. The second consequence, noisy statistics at small batch sizes, is the failure mode covered in the last section.

With the statistics frozen, BN at inference is a fixed affine map per feature, and Algorithm 2 folds it together with $\gamma$ and $\beta$:

The folded form costs one multiply and one add per value, and it can be merged into the weights and bias of the preceding linear or convolutional layer, so a deployed BN layer costs nothing extra. The paper notes that moving averages kept during training can replace the after-training average in (6), which lets you measure validation accuracy while the model trains; libraries do this. PyTorch keeps a running average with momentum 0.1 and feeds it the unbiased per-batch variance, applying $m/(m-1)$ batch by batch; Keras keeps one with momentum 0.99 of the biased variance and skips the correction. The listing below follows PyTorch and uses (5) for the backward pass.

# Batch norm for ONE feature over a mini-batch x of shape (m,), with m > 1.
# (For a conv layer, x holds one channel's values over batch and positions.)
def bn_forward(x, gamma, beta, run_mean, run_var, training,
               momentum=0.1, eps=1e-5):
    if training:
        m = x.shape[0]
        mu = x.mean()                          # (1/m) sum x_i
        var = ((x - mu) ** 2).mean()           # biased: divide by m
        run_mean = (1 - momentum) * run_mean + momentum * mu
        run_var = (1 - momentum) * run_var + momentum * var * m / (m - 1)
    else:
        mu, var = run_mean, run_var            # frozen population stats
    s = (var + eps) ** 0.5
    x_hat = (x - mu) / s                       # mean 0, variance ~1
    return gamma * x_hat + beta, x_hat, s, run_mean, run_var

def bn_backward(dy, gamma, x_hat, s):
    # dy = dloss/dy, shape (m,). Equation (5); eps is inside s.
    g = dy * gamma                             # dloss/dx_hat
    dx = (g - g.mean() - x_hat * (g * x_hat).mean()) / s
    return dx, (dy * x_hat).sum(), dy.sum()    # dx, dgamma, dbeta

Why batch norm allows higher learning rates

In the paper, BN networks train at learning rates 5 and 30 times higher than the original network, while a 5× rate on the original Inception without BN drove its parameters to "machine infinity". Section 3.3 gives an argument from scale. Multiply a layer's weights by a constant $a > 0$ and every $x_i$ in the batch is multiplied by $a$, so the mean and standard deviation scale by $a$ too and the normalized output does not change:

The gradient passed to the layer below does not depend on the weight scale, and the gradient on the weights shrinks as $1/a$: doubling the weights halves their gradient. The paper concludes that "larger weights lead to smaller gradients, and Batch Normalization will stabilize the parameter growth".

Two further facts follow from scale invariance. Because the loss does not change when $W$ is scaled, its gradient has no component along $W$ and is perpendicular to it. A gradient step perpendicular to $W$ can only lengthen it; by Pythagoras, with $\nabla$ the gradient and $\eta$ the learning rate:

A step of length $\eta\|\nabla\|$ perpendicular to a vector of length $\|W\|$ turns it by an angle of about $\eta\|\nabla\|/\|W\|$, where $\|\nabla\|$ itself is proportional to $1/\|W\|$. The direction, the only thing the loss depends on, changes at an effective rate of $\eta/\|W\|^2$. A learning rate that is too large makes a big step, which lengthens $W$, which lowers the effective rate for every step after. Later work develops this reading: van Laarhoven (2017) showed that L2 weight decay has no regularizing effect on normalized weights and acts through their scale on the effective learning rate, and Arora, Li and Lyu (2018) analyze BN as automatic tuning of the learning rate.

The figure shows a two-weight unit whose loss depends only on the direction of $w$, with the best direction dashed. At $\eta = 1.5$ from $\|w\| = 0.5$, the first step overshoots, lengthens $w$ to 1.63, and cuts the effective step from 6.0 to 0.57, after which the angle converges. Push $\eta$ to 8: one step takes $\|w\|$ from 0.5 to 8.26, the effective step drops to 0.117, and the descent still converges, slowly. Then try a small $\eta$ with a large starting $\|w\|$, where it barely moves.

Section 3.3 also offers a second, explicitly hedged argument: "We further conjecture that Batch Normalization may lead the layer Jacobians to have singular values close to 1", which would keep gradient magnitudes steady through depth. The derivation assumes the normalized inputs are Gaussian and uncorrelated and the layer is linear, and the paper says the precise effect "remains an area of further study". Yang et al. (2019) settled it the other way for networks at initialization: in the infinite-width limit, BN pushes the singular values away from 1 and gradients explode with depth.

A third view, from Santurkar et al. (2018) and covered in the last section, is that BN makes the loss surface smoother. For intuition on why smoothness allows a larger step: gradient descent on a quadratic with curvature $\kappa$ is stable only while the step size is below $2/\kappa$. The two toy one-parameter landscapes below share a learning rate; the left one has high curvature, the right one low.

At the default rate of 0.45 the left ball bounces between walls and the right one settles. At the maximum, 1.40, both diverge: smoothing raises the largest usable learning rate, and some rate is still too large.

ImageNet results: 14× fewer training steps

The testbed is a variant of Inception, the network that won the 2014 ImageNet classification challenge, on the 1000-class ImageNet task. The variant has 13.6 million parameters, replaces each 5×5 convolution with two 3×3 convolutions, uses ReLU, and has no fully connected layers other than the final softmax. It is trained by SGD with momentum at batch size 32 with an initial learning rate of 0.0015. BN goes before every nonlinearity, applied convolutionally, and the rest of the architecture is unchanged, so the variants differ only in BN and training settings.

Adding BN and nothing else (BN-Baseline) reaches Inception's final 72.2% top-1 validation accuracy (the fraction of held-out images whose highest-scoring class is correct, from a single crop) in 13.3 million steps instead of 31.0 million and tops out at 72.7%. The authors then changed the training setup to use what BN allows. The initial learning rate went up 5×, to 0.0075. Dropout was removed. The L2 weight penalty was reduced by a factor of 5. The learning rate was decayed 6 times faster. Local Response Normalization, an older cross-channel normalization in Inception, was removed. Training examples were shuffled more thoroughly, and photometric distortions of the training images were reduced. That model, BN-x5, reaches 72.2% in 2.1 million steps, 7% of Inception's 31.0 million: the 14× in the abstract. Its best accuracy is 73.0%. The comparison is in training steps; the paper does not report wall-clock time or the extra cost of a BN step.

BN-x30 uses an initial rate of 0.045, 30× Inception's. It trains somewhat slower than BN-x5 at first, needing 2.7 million steps to reach 72.2%, but it keeps improving and reaches 74.8% after 6 million steps. Drag the target line in the figure: at 72.2% the crossing points reproduce the step counts, and above 73.0% only BN-x30 still crosses.

The sigmoid experiment tests the saturation argument directly. BN-x5-Sigmoid, the BN-x5 network with sigmoid in place of ReLU, reaches 69.8%. The original Inception with sigmoid and no BN "never achieves better than 1/1000 accuracy", which is chance for 1000 classes.

An ensemble of 6 networks, each based on BN-x30 with small modifications (among them Dropout at 5% or 10%), reaches 4.9% top-5 validation error and 4.82% top-5 error on the 100,000-image test set, as scored by the ILSVRC server. Top-5 error is the fraction of images whose correct label is not among the model's 5 highest-scoring classes. The previous results in the paper's comparison table are the GoogLeNet ensemble at 6.67% and Deep Image at 5.98%; the text also cites He et al. (2015) at 4.94%. The paper says the ensemble exceeds "the estimated accuracy of human raters", citing Russakovsky et al. (2014). That estimate is one trained annotator's 5.1% top-5 error on 1,500 test images; a second annotator with less training scored about 12.0% on 258 images, and Russakovsky et al. expected humans to stay ahead of the models of the time "only by use of significant effort, expertise, and time".

What later work found about why batch norm works

The paper's title states its explanation: BN accelerates training by reducing internal covariate shift. In 2018, Santurkar, Tsipras, Ilyas and Madry tested that explanation in three ways.

First, they trained BN networks with random noise added after every BN layer, drawn from a distribution with non-zero mean and non-unit variance that changed at every step. The resulting activations were less stable than those of a network without BN at all, and the noisy-BN network trained almost as well as the ordinary BN network, while the same noise added to a network without BN stopped it from training entirely. Second, they defined internal covariate shift through how much a layer's gradient changes when the layers before it are updated, and found that BN networks did not show less of it, and sometimes more; they call the link between BN's performance and internal covariate shift "tenuous, at best". Third, they measured what BN does change: the loss and its gradient vary more slowly with the parameters (better Lipschitz constants and smoothness), so a gradient computed at one point still points in a useful direction over a longer step. They also found that normalizing by the $\ell_1$, $\ell_2$ or $\ell_\infty$ norm, none of which produces mean-0, variance-1 activations, trains comparably to BN. Their smoothness theorems hold for specific simplified settings, and their result where $\ell_1$ normalization beats BN is for deep linear networks.

Santurkar et al. do not claim the gains are illusory or that activation distributions never drift. Their experiments show that keeping the distributions stable is not what makes BN networks train faster, and their smoothness results, together with the scale-invariance argument of Figure 6, account for the higher usable learning rates.

The paper does not discuss batch size. With small batches the statistics are noisy, as Figure 5 showed, and the difference between batch statistics in training and population statistics at inference grows. Wu and He (2018) measured ResNet-50 on ImageNet: validation error with BN rose from 23.6% at 32 images per GPU to 34.7% at 2 images per GPU, while their Group Normalization, which normalizes over groups of channels within each example, stayed at 24.1%. Ioffe's own Batch Renormalization (2017) corrects the batch statistics toward the moving averages during training. Layer Normalization (Ba, Kiros and Hinton, 2016) normalizes over all features of one example, with no batch involved; the Transformer adopted it the following year, and Transformer models have used it since. In convolutional networks BN stayed standard: ResNet applies it after every convolution.