NumPy Broadcasting in Python: Rules and Common Bugs (Part 3)

NumPy broadcasting is the set of rules that lets arithmetic work on arrays of different shapes, so you can write arr + 10 instead of a loop. NumPy compares shapes from the right, and each pair of dimensions must be equal or contain a 1. Those two rules cover the whole system.

The danger is that broadcasting can succeed when you wanted it to fail. Add a (3,) array to a (3, 1) array and you get a 3x3 matrix, not an error. This part of the series is for readers who know indexing and slicing from Part 2. It works through the rules with shape diagrams, the silent bug, the tools that control alignment, and two worked examples. Outputs come from NumPy 2.2.

What are the NumPy broadcasting rules?

When NumPy operates on two arrays, it compares their shapes element by element, starting from the trailing (rightmost) dimension. According to the official broadcasting guide, two dimensions are compatible when:

1. they are equal, or
2. one of them is 1

If a dimension is missing from the shorter shape, it is treated as 1. If any pair fails both tests, then the operation raises ValueError: operands could not be broadcast together.

Worked through:

A      (3, 4)
B         (4,)
-----------------
       (3, 4)   ✓  B is treated as (1,4), stretched down 3 rows

A      (3, 1)
B         (4,)
-----------------
       (3, 4)   ✓  A stretched across, B stretched down

A      (3, 4)
B         (3,)
-----------------
        error   ✗  4 vs 3 in the trailing dimension

You can also ask NumPy to apply the rules without doing any arithmetic:

import numpy as np

np.broadcast_shapes((3, 1), (4,))            # (3, 4)
np.broadcast_shapes((8, 1, 6, 1), (7, 1, 5))  # (8, 7, 6, 5)

How does broadcasting work in practice?

In practice, a scalar, a row or a column is stretched across a matrix so one operation applies everywhere. Scalar with array is the simplest broadcast:

a = np.arange(12).reshape(3, 4)
a + 100        # every element gets 100 added
a * 2

Next, a row vector across a matrix:

matrix = np.array([[1, 2, 3],
                   [4, 5, 6]])          # (2, 3)
row = np.array([10, 20, 30])            # (3,)

matrix + row
# [[11 22 33]
#  [14 25 36]]

The row is applied to each row of the matrix. This is the pattern behind per-column operations such as centring:

data = np.array([[1., 2., 3.],
                 [4., 5., 6.],
                 [7., 8., 9.]])

col_means = data.mean(axis=0)     # [4. 5. 6.], shape (3,)
centred = data - col_means        # each column gets its own mean subtracted

A column vector down a matrix, by contrast, needs an explicit second dimension:

col = np.array([[10],
                [20]])            # (2, 1)

matrix + col
# [[11 12 13]
#  [24 25 26]]

Why does NumPy broadcasting fail silently?

Because shapes can be compatible in a way you did not intend. An (n,) array and an (n, 1) array both stretch, so NumPy returns an (n, n) matrix instead of raising an error, and nothing warns you.

a = np.array([1, 2, 3])            # shape (3,)
b = np.array([[10], [20], [30]])   # shape (3, 1)

a + b
# [[11 12 13]
#  [21 22 23]
#  [31 32 33]]

You expected [11 22 33]. Instead, you got a 3x3 matrix.

Right-aligning the shapes, (3,) becomes (1, 3) and meets (3, 1). Both dimensions have a 1, so both stretch, and the result is the outer sum.

Nothing errors. Downstream you get wrong numbers, or a memory explosion when the arrays are large: a (10000, 1) with a (10000,) produces a hundred-million-element array.

The fix is always to check shapes:

print(a.shape, b.shape)   # (3,) (3, 1)
a + b.ravel()             # [11 22 33]  - flatten b to (3,)

The row-wise version of this bug raises an error instead, which is the lucky case:

M = np.array([[1., 2., 3.],
              [4., 5., 6.]])

M.mean(axis=1)                          # [2. 5.], shape (2,)
M - M.mean(axis=1, keepdims=True)       # [[-1. 0. 1.] [-1. 0. 1.]]

Without keepdims=True, M - M.mean(axis=1) fails with shapes (2,3) and (2,), because the trailing dimensions 3 and 2 do not match.

How do you control broadcasting with np.newaxis?

Add a length-one dimension exactly where you want the stretch to happen, with np.newaxis, reshape or keepdims=True. The NumPy indexing guide describes it as a way to insert a new axis of length one, and it is simply an alias for None:

a = np.array([1, 2, 3])

a[:, np.newaxis].shape      # (3, 1) - column
a[np.newaxis, :].shape      # (1, 3) - row
a.reshape(-1, 1).shape      # (3, 1) - same as the first

This makes an outer operation explicit rather than accidental:

x = np.array([1, 2, 3])
y = np.array([10, 20])

x[np.newaxis, :] * y[:, np.newaxis]
# [[10 20 30]
#  [20 40 60]]

Written that way, the intent is visible to whoever reads it next.

Worked example: normalising features

Broadcasting makes standard preprocessing a two-line operation:

rng = np.random.default_rng(0)
X = rng.normal(size=(100, 5))

means = X.mean(axis=0)      # (5,)
stds  = X.std(axis=0)       # (5,)

X_scaled = (X - means) / stds
X_scaled.std(axis=0).round(10)   # [1. 1. 1. 1. 1.]

Each of the 5 columns gets its own mean and standard deviation, applied down all 100 rows, with no loop and no copy of the statistics.

The axis argument is the other half of this. axis=0 collapses rows, giving one value per column, while axis=1 collapses columns, giving one value per row. Part 4 covers that argument in detail.

Worked example: pairwise distances without a loop

Adding a dimension to each copy of a point array turns every pair into one subtraction:

pts = np.array([[0, 0], [3, 4], [6, 8]])       # (3, 2)

diff = pts[:, np.newaxis, :] - pts[np.newaxis, :, :]   # (3, 3, 2)
dist = np.sqrt((diff ** 2).sum(axis=-1))
# [[ 0.  5. 10.]
#  [ 5.  0.  5.]
#  [10.  5.  0.]]

The trade-off is memory. For n points the intermediate array is n × n × 2, so the NumPy documentation itself notes that broadcasting can be a bad idea when it creates a very large intermediate array. For tens of thousands of points, process in chunks instead.

NumPy broadcasting quick reference

Common NumPy broadcasting mistakes

• Mixing (n,) and (n, 1). The sum becomes an (n, n) matrix; check .shape first.
• Dropping keepdims=True on row-wise statistics. The result no longer lines up with the rows.
• Broadcasting in place into a smaller array. z = np.zeros(3); z += np.ones((3, 3)) raises "non-broadcastable output operand", because the result cannot fit in z.
• Building huge intermediates. Pairwise tricks on large inputs can use more memory than a chunked loop.

Related reading

NumPy Part 2: indexing and slicing precedes this, and NumPy Part 1 starts the series. NumPy Part 4: aggregation and linear algebra continues. For the operators themselves, see Python operators.