Linear Algebra for Machine Learning

The four ideas that make up linear algebra for machine learning

Linear algebra for machine learning is a much smaller subject than a university course suggests. You need vectors, matrices, the dot product and matrix multiplication. Add norms and a rough sense of what eigenvectors are and you can read most papers in the field.

You do not need to invert matrices by hand, compute determinants, or prove anything about vector spaces. Those skills are for people who write the libraries. Your job is to know what shape your data is, what operation the model is performing on it, and why an error message about (400, 3) and (4,) means what it means.

Everything here uses numpy. A vector is a one-dimensional array, a matrix is two-dimensional, and the same operators work on both.

The dot product is the model

Take the delivery-time formula from what machine learning actually is: about three minutes per kilometre, one and three-quarter minutes per kilogram, nine minutes of fixed overhead. Predicting one delivery means multiplying each feature by its weight and adding the results.

That operation is the dot product. Two vectors of the same length, multiplied element by element, then summed:

import numpy as np

weights = np.array([3.04, 1.73])    # minutes per km, minutes per kg
parcel  = np.array([18.0, 2.5])     # an 18 km delivery weighing 2.5 kg
bias    = 8.76

print(round(float(np.dot(weights, parcel) + bias), 2))   # 67.81

Linear regression, logistic regression, support vector machines and every neuron in a neural network are doing this. The differences lie in how the weights are found and what happens to the result afterwards. The core arithmetic never changes: multiply features by weights, add them up.

Write the dot product as a formula and it is w₁x₁ + w₂x₂ + ... + wₙxₙ — each weight times its matching feature, summed across all n features. That is the entire mathematical content of a linear model.

Matrix multiplication does it for every row at once

One prediction is a dot product. A thousand predictions is a matrix multiplication, and this is where the speed comes from.

X = np.array([[2.0, 0.5],      # each row: one delivery
              [5.0, 1.2],
              [9.0, 0.8],
              [12.0, 3.5]])

print(X.shape, weights.shape)        # (4, 2) (2,)
print((X @ weights + bias).round(2)) # [15.7  26.04 37.5  51.29]

The @ operator is matrix multiplication. One expression, four predictions, and it would be four million predictions with no change to the code — the loop runs inside compiled BLAS routines rather than in Python.

Shapes are the thing that breaks

Almost every numpy error in machine learning is a shape error, and the rule is worth committing to memory. To multiply A @ B, the last dimension of A must equal the first dimension of B. (4, 2) @ (2,) works and gives (4,). (4, 2) @ (4,) raises ValueError.

The convention in scikit-learn is that X is always (n_samples, n_features) — rows are observations, columns are variables, which is the same X and y contract described in supervised learning. Get that transposed and either nothing works or, worse, it works and trains on nonsense. When a shape error appears, print .shape on both operands before changing anything. A.T gives the transpose, flipping rows and columns.

Solving for the weights directly

You can also go the other way: given X and the answers y, find the weights. For linear regression there is a closed-form solution, the normal equation, and it is one line of numpy.

y  = np.array([14., 26., 38., 52.])
Xb = np.hstack([np.ones((4, 1)), X])          # column of 1s for the intercept

theta = np.linalg.solve(Xb.T @ Xb, Xb.T @ y)
print(theta.round(4))    # [6.9984 3.3164 1.5246]

from sklearn.linear_model import LinearRegression
model = LinearRegression().fit(X, y)
print(round(model.intercept_, 4), model.coef_.round(4))
# 6.9984 [3.3164 1.5246]

Identical to four decimal places, because LinearRegression solves the same system. Use np.linalg.solve rather than computing an inverse with np.linalg.inv — it is faster and numerically better behaved. If Xb.T @ Xb turns out to be singular, which happens when two columns are perfectly correlated, solve raises an error and np.linalg.lstsq gives you the least-squares answer instead.

This route stops being practical around tens of thousands of features, because the cost grows roughly with the cube of the feature count. That is the point at which iterative methods take over — the guess, measure, adjust loop behind model training — which is what the rest of this section builds towards.

Norms measure the size of a vector

A norm collapses a vector to a single non-negative number describing how big it is. Two matter.

The L2 norm is ordinary straight-line length: square each element, add, take the square root. The L1 norm is the sum of absolute values.

w_spread = np.array([0.9, 0.8, 0.85, 0.7])   # four moderate weights
w_spike  = np.array([1.6, 0.0, 0.0, 0.0])    # one large, three zero

for w in (w_spread, w_spike):
    print(round(float(np.linalg.norm(w, 1)), 3),
          round(float(np.linalg.norm(w, 2)), 3))
# 3.25 1.632
# 1.6  1.6

Read that carefully, because it explains a great deal about regularisation. Under the L1 norm the spread-out weights score 3.25 and the spiky one 1.6, so a penalty based on L1 prefers the sparse solution and pushes coefficients to exactly zero. Under L2 they score 1.632 and 1.6, almost the same, so an L2 penalty shrinks everything a little and zeroes nothing. Ridge and lasso regression are that difference and nothing more.

Eigenvectors are the one item on that list you can defer. They matter for principal component analysis, and the intuition — directions along which a matrix only stretches, never rotates — is enough until you reach that chapter.

The full reference for these operations is numpy’s linear algebra routines, and the solver used above is documented under numpy.linalg.solve.

Frequently Asked Questions

How much linear algebra do you need for machine learning?

Less than a full course. Vectors, matrices, the dot product, matrix multiplication and shape rules cover the day-to-day work. Norms help you understand regularisation, and eigenvectors matter for dimensionality reduction. Determinants, proofs and hand computation of inverses are not needed for applied work.

Why does matrix multiplication matter so much in machine learning?

It applies the same weights to every row of your data in one operation, executed by optimised compiled libraries rather than a Python loop. That single fact is why a model can score a million records in under a second. Nearly every model’s forward pass is a matrix product.

What is the difference between the L1 and L2 norm?

L1 sums absolute values; L2 takes the square root of summed squares. L1 treats a single large weight as cheaper than several moderate ones, so it drives coefficients to exactly zero and produces sparse models. L2 penalises large weights more steeply and shrinks all of them without eliminating any.

Key Takeaways

• Print .shape on both operands whenever numpy raises a dimension error, because almost every linear algebra bug in machine learning is a shape mismatch rather than a mathematical mistake.
• Keep X as (n_samples, n_features) throughout your code, since a transposed design matrix can train silently on nonsense.
• Reach for @ instead of writing loops over rows, as the vectorised form is both shorter and orders of magnitude faster.
• Use np.linalg.solve rather than inverting a matrix explicitly, and fall back to np.linalg.lstsq when columns are collinear.
• Remember that L1 produces sparsity and L2 produces shrinkage, because that one distinction explains the practical difference between lasso and ridge.