Linear Algebra for Machine Learning: Vectors, Matrices, and Transformations

Linear algebra gives machine learning a compact way to represent data, express models, and perform many calculations at once. A photograph becomes a grid of numbers, a sentence becomes a sequence of embedding vectors, and an entire dataset becomes a matrix whose rows describe examples. Training then changes vectors of parameters so that matrix-based predictions better match observed targets. This guide develops the essential ideas from geometry and arithmetic, connects them to real AI systems, and emphasizes the shape reasoning that makes implementations reliable.

Why linear algebra is the language of AI

Machine learning begins by turning observations into numbers. Suppose a house is represented by floor area, age, and distance to a station. One house can be encoded as a vector such as

$$\mathbf{x}= \begin{bmatrix} 120\\ 8\\ 0.7 \end{bmatrix}.$$

A dataset of one thousand houses is naturally represented by a matrix with one thousand rows and three columns. A linear prediction model stores one weight per feature in another vector. Instead of writing a separate equation for every house, we multiply the data matrix by the weight vector and obtain every prediction simultaneously.

That compression is not merely convenient notation. Modern processors, graphics processing units, and specialized AI accelerators are built to execute large blocks of matrix arithmetic efficiently. The mathematical abstraction and the hardware implementation reinforce one another: a model can be described in a few equations, then evaluated in parallel across millions of values.

Linear algebra appears throughout the machine-learning pipeline:

• Data representation: tables, images, audio windows, token embeddings, and batches are arrays with meaningful dimensions.
• Model evaluation: linear regression, dense neural-network layers, convolutions, and attention all rely heavily on multiplication and addition over arrays.
• Similarity: recommendation and retrieval systems compare embedding vectors with dot products or cosine similarity.
• Dimensionality reduction: methods such as principal component analysis identify important directions in a dataset.
• Optimization: gradients have the same shapes as the parameters they update, while curvature can be represented by matrices.
• Statistics: covariance matrices describe how features vary together.

Even nonlinear models are usually assembled from linear operations and simple nonlinear functions. A neural-network layer computes an affine map

$$\mathbf{z}=\mathbf{W}\mathbf{x}+\mathbf{b}$$

and then applies a nonlinear activation such as ReLU:

$$\operatorname{ReLU}(z_i)=\max(0,z_i).$$

Without the nonlinear activation, stacking many layers would still collapse into one linear map. With it, neural networks can model curved and piecewise relationships, but the expensive core remains dominated by matrix operations.

Learning linear algebra therefore serves two goals. It explains what a model is doing conceptually, and it helps you reason about code. Statements such as “this tensor has shape batch by features,” “these embeddings are nearly orthogonal,” or “the design matrix is rank deficient” become precise rather than mysterious.

Vectors, norms, and geometric intuition

What a vector represents

A vector is an ordered list of numbers. Depending on context, it can represent a point, a direction, a displacement, a feature record, or a set of model parameters. We usually write

$$\mathbf{x}= \begin{bmatrix} x_1\\ x_2\\ \vdots\\ x_n \end{bmatrix} \in \mathbb{R}^n.$$

The notation $\mathbb{R}^n$ means that $\mathbf{x}$ has $n$ real-valued components. Order matters: if the coordinates represent height, weight, and age, exchanging the first and third entries changes the meaning.

In two or three dimensions, a vector can be drawn as an arrow from the origin. Its coordinates tell us how far the arrow travels along each axis. A vector $\begin{bmatrix}3 & 2\end{bmatrix}^T$ moves three units horizontally and two vertically. This picture remains useful in higher dimensions even though we cannot draw them directly.

Vector addition combines displacements component by component:

$$\begin{bmatrix}2\\1\end{bmatrix} + \begin{bmatrix}-1\\3\end{bmatrix} = \begin{bmatrix}1\\4\end{bmatrix}.$$

Multiplication by a scalar stretches or reverses a vector:

$$-2 \begin{bmatrix}2\\1\end{bmatrix} = \begin{bmatrix}-4\\-2\end{bmatrix}.$$

These two operations—addition and scalar multiplication—are the foundation of linear combinations. Given vectors $\mathbf{v}_1,\ldots,\mathbf{v}_k$, any expression

$$c_1\mathbf{v}_1+\cdots+c_k\mathbf{v}_k$$

is a linear combination of those vectors.

Norms measure size

A norm assigns a nonnegative size to a vector. The Euclidean or $L_2$ norm is the familiar straight-line length:

$$\|\mathbf{x}\|_2=\sqrt{x_1^2+x_2^2+\cdots+x_n^2}.$$

For $\mathbf{x}=\begin{bmatrix}3 & 4\end{bmatrix}^T$, the $L_2$ norm is $5$. Distances follow directly because the distance from $\mathbf{x}$ to $\mathbf{y}$ is $\|\mathbf{x}-\mathbf{y}\|_2$.

Other norms emphasize different properties:

A unit vector has norm one. Normalizing a nonzero vector preserves its direction while changing its length to one:

$$\hat{\mathbf{x}}=\frac{\mathbf{x}}{\|\mathbf{x}\|_2}.$$

Normalization is useful when direction matters more than scale. For example, two document embeddings may point in similar directions even if one has a larger magnitude. A warning is essential: the zero vector cannot be normalized because division by $\|\mathbf{0}\|=0$ is undefined. Production code commonly adds a small tolerance or handles zero vectors explicitly.

Norms also appear in training objectives. Ridge regression penalizes $\|\mathbf{w}\|_2^2$, discouraging very large weights smoothly. Lasso penalizes $\|\mathbf{w}\|_1$, which can drive some weights exactly to zero. The choice of norm therefore changes not only measurement but also learned model behavior.

Matrices as linear maps

From rectangular arrays to transformations

A matrix is a rectangular array of numbers with $m$ rows and $n$ columns:

$$\mathbf{A}\in\mathbb{R}^{m\times n}.$$

It can store data, but a deeper interpretation is that $\mathbf{A}$ defines a function from $\mathbb{R}^n$ to $\mathbb{R}^m$:

$$\mathbf{y}=\mathbf{A}\mathbf{x}.$$

The input must have $n$ components, and the output has $m$ components. If

$$\mathbf{A}= \begin{bmatrix} a_{11} & a_{12}\\ a_{21} & a_{22} \end{bmatrix}, \qquad \mathbf{x}= \begin{bmatrix} x_1\\ x_2 \end{bmatrix},$$

then

$$\mathbf{A}\mathbf{x}= \begin{bmatrix} a_{11}x_1+a_{12}x_2\\ a_{21}x_1+a_{22}x_2 \end{bmatrix}.$$

The transformation is called linear because it preserves addition and scalar multiplication:

$$\mathbf{A}(\mathbf{x}+\mathbf{y})=\mathbf{A}\mathbf{x}+\mathbf{A}\mathbf{y},$$ $$\mathbf{A}(c\mathbf{x})=c\mathbf{A}\mathbf{x}.$$

Linear maps can rotate, reflect, scale, shear, project, or combine these effects. They always map the origin to the origin. A transformation $\mathbf{A}\mathbf{x}+\mathbf{b}$ with nonzero bias $\mathbf{b}$ is technically affine rather than linear, because it shifts the origin. Machine-learning discussions often use “linear layer” for an affine layer, so remembering this distinction prevents conceptual confusion.

Another useful interpretation comes from the columns of $\mathbf{A}$. Write

$$\mathbf{A}= \begin{bmatrix} \vert & \vert & & \vert\\ \mathbf{a}_1 & \mathbf{a}_2 & \cdots & \mathbf{a}_n\\ \vert & \vert & & \vert \end{bmatrix}.$$

Then

$$\mathbf{A}\mathbf{x}=x_1\mathbf{a}_1+x_2\mathbf{a}_2+\cdots+x_n\mathbf{a}_n.$$

The output is a linear combination of the matrix columns, with the input coordinates supplying the coefficients. This column view will make span, rank, and solvability much easier to understand.

Composition and order

Suppose $\mathbf{B}$ transforms $\mathbf{x}$ first and $\mathbf{A}$ transforms the result:

$$\mathbf{x}\xrightarrow{\mathbf{B}}\mathbf{B}\mathbf{x} \xrightarrow{\mathbf{A}}\mathbf{A}(\mathbf{B}\mathbf{x}).$$

The combined transformation is

$$\mathbf{A}\mathbf{B}\mathbf{x}.$$

Matrix multiplication therefore represents composition. The rightmost operation occurs first, just as in ordinary function composition. This order matters because matrix multiplication is generally not commutative:

$$\mathbf{A}\mathbf{B}\ne\mathbf{B}\mathbf{A}.$$

Rotating and then stretching along the horizontal axis usually produces a different result from stretching first and then rotating. In neural networks, changing layer order likewise changes the function.

The dimensions encode whether composition is possible. If $\mathbf{B}\in\mathbb{R}^{p\times n}$ and $\mathbf{A}\in\mathbb{R}^{m\times p}$, then $\mathbf{A}\mathbf{B}\in\mathbb{R}^{m\times n}$. The shared inner dimension $p$ must agree because $\mathbf{B}$ produces $p$ values and $\mathbf{A}$ expects $p$ inputs.

Matrix multiplication and neural networks

The row-by-column rule

For $\mathbf{A}\in\mathbb{R}^{m\times n}$ and $\mathbf{B}\in\mathbb{R}^{n\times p}$, the product $\mathbf{C}=\mathbf{A}\mathbf{B}$ has shape $m\times p$. Its entries are

$$c_{ij}=\sum_{k=1}^{n}a_{ik}b_{kj}.$$

Each output entry is the dot product between row $i$ of $\mathbf{A}$ and column $j$ of $\mathbf{B}$. A practical shape mnemonic is

$$(m\times \color{#555}{n})(\color{#555}{n}\times p)\longrightarrow(m\times p).$$

The inner dimensions must match; the outer dimensions remain.

Matrix multiplication can also be understood as processing several vectors at once. The columns of $\mathbf{B}$ are separate input vectors. Multiplying by $\mathbf{A}$ applies the same transformation to each column:

$$\mathbf{A} \begin{bmatrix} \vert & \vert\\ \mathbf{b}_1 & \mathbf{b}_2\\ \vert & \vert \end{bmatrix} = \begin{bmatrix} \vert & \vert\\ \mathbf{A}\mathbf{b}_1 & \mathbf{A}\mathbf{b}_2\\ \vert & \vert \end{bmatrix}.$$

Why dense layers are matrix multiplications

Consider a dense neural-network layer with $d_{\text{in}}$ input features and $d_{\text{out}}$ output units. For one column-vector input,

$$\mathbf{z}=\mathbf{W}\mathbf{x}+\mathbf{b},$$

where

$$\mathbf{W}\in\mathbb{R}^{d_{\text{out}}\times d_{\text{in}}}, \quad \mathbf{x}\in\mathbb{R}^{d_{\text{in}}}, \quad \mathbf{b}\in\mathbb{R}^{d_{\text{out}}}.$$

Every output neuron computes a weighted sum of all input features. Each row of $\mathbf{W}$ stores one neuron's weights, so evaluating all neurons is exactly a matrix-vector multiplication.

Libraries often organize a batch as rows rather than columns. If $\mathbf{X}$ has shape $(B,d_{\text{in}})$ for batch size $B$, the equivalent expression is

$$\mathbf{Z}=\mathbf{X}\mathbf{W}^T+\mathbf{b},$$

with output shape $(B,d_{\text{out}})$. The bias is broadcast across all $B$ rows.

Transformers are also filled with matrix multiplication. Queries, keys, and values are obtained by learned projections such as $\mathbf{Q}=\mathbf{X}\mathbf{W}_Q$. Attention scores are built from

$$\frac{\mathbf{Q}\mathbf{K}^T}{\sqrt{d_k}},$$

and the resulting weights multiply $\mathbf{V}$. Convolutions can similarly be implemented as structured matrix multiplications, even though libraries use more specialized kernels. Calling modern AI “mostly matmuls” is an oversimplification, but it captures the dominant arithmetic in many models.

Transpose, identity, and inverse

Transpose

The transpose flips rows and columns:

$$\mathbf{A}= \begin{bmatrix} 1&2&3\\ 4&5&6 \end{bmatrix} \quad\Longrightarrow\quad \mathbf{A}^T= \begin{bmatrix} 1&4\\ 2&5\\ 3&6 \end{bmatrix}.$$

If $\mathbf{A}$ is $m\times n$, then $\mathbf{A}^T$ is $n\times m$. Transposition obeys several useful rules:

$$(\mathbf{A}^T)^T=\mathbf{A}, \qquad (\mathbf{A}+\mathbf{B})^T=\mathbf{A}^T+\mathbf{B}^T,$$ $$(\mathbf{A}\mathbf{B})^T=\mathbf{B}^T\mathbf{A}^T.$$

Notice the reversed order in the product rule. That reversal follows from swapping the role of rows and columns.

A square matrix is symmetric when $\mathbf{A}=\mathbf{A}^T$. Covariance matrices and matrices such as $\mathbf{X}^T\mathbf{X}$ are symmetric. Symmetry brings valuable properties, including real eigenvalues and orthogonal eigenvectors under common conditions.

Identity and inverse

The identity matrix $\mathbf{I}$ has ones on its main diagonal and zeros elsewhere. It leaves vectors unchanged:

$$\mathbf{I}\mathbf{x}=\mathbf{x}, \qquad \mathbf{A}\mathbf{I}=\mathbf{I}\mathbf{A}=\mathbf{A}.$$

An inverse of a square matrix $\mathbf{A}$ is a matrix $\mathbf{A}^{-1}$ satisfying

$$\mathbf{A}^{-1}\mathbf{A} = \mathbf{A}\mathbf{A}^{-1} =\mathbf{I}.$$

If $\mathbf{y}=\mathbf{A}\mathbf{x}$ and an inverse exists, then

$$\mathbf{x}=\mathbf{A}^{-1}\mathbf{y}.$$

Geometrically, an invertible transformation loses no information: every output corresponds to exactly one input. A rotation by $30^\circ$ is invertible because rotating by $-30^\circ$ reverses it. Scaling one axis by two is invertible because scaling that axis by one half reverses it.

When inverses fail

Not every matrix has an inverse. A non-square matrix cannot have a standard two-sided inverse. Even a square matrix is singular if it compresses different inputs to the same output. For example,

$$\mathbf{A}= \begin{bmatrix} 1&2\\ 2&4 \end{bmatrix}$$

maps the plane onto a line. Its second row is twice the first, so one equation contains no new information. No transformation can reconstruct the lost dimension.

For a $2\times2$ matrix,

$$\mathbf{A}= \begin{bmatrix} a&b\\ c&d \end{bmatrix},$$

an inverse exists exactly when the determinant $ad-bc$ is nonzero. When it exists,

$$\mathbf{A}^{-1}= \frac{1}{ad-bc} \begin{bmatrix} d&-b\\ -c&a \end{bmatrix}.$$

An inverse can also be numerically dangerous when the determinant is nonzero but the matrix is nearly singular. Small changes in data may then produce large changes in the computed solution. This behavior is summarized by the condition number: a large condition number indicates sensitivity.

In numerical machine learning, explicitly computing $\mathbf{A}^{-1}$ is rarely the best way to solve $\mathbf{A}\mathbf{x}=\mathbf{b}$. A solver based on LU, QR, Cholesky, or singular value decomposition is typically faster and more stable. Use numpy.linalg.solve(A, b) for a square system instead of np.linalg.inv(A) @ b. For rectangular or least-squares problems, use numpy.linalg.lstsq.

Rank, linear independence, and span

Span describes reachable directions

The span of vectors $\mathbf{v}_1,\ldots,\mathbf{v}_k$ is the set of all their linear combinations:

$$\operatorname{span}(\mathbf{v}_1,\ldots,\mathbf{v}_k) = \left\{ c_1\mathbf{v}_1+\cdots+c_k\mathbf{v}_k \,:\, c_i\in\mathbb{R} \right\}.$$

One nonzero vector in $\mathbb{R}^2$ spans a line through the origin. Two nonparallel vectors span the whole plane. In $\mathbb{R}^3$, two nonparallel vectors usually span a plane through the origin, while three suitably independent vectors span all three dimensions.

Because $\mathbf{A}\mathbf{x}$ is a linear combination of the columns of $\mathbf{A}$, the set of all possible outputs is the column space:

$$\operatorname{Col}(\mathbf{A}) = \operatorname{span}(\mathbf{a}_1,\ldots,\mathbf{a}_n).$$

The equation $\mathbf{A}\mathbf{x}=\mathbf{b}$ has an exact solution only if $\mathbf{b}$ lies in this column space.

Independence detects redundancy

Vectors are linearly independent if the equation

$$c_1\mathbf{v}_1+\cdots+c_k\mathbf{v}_k=\mathbf{0}$$

has only the trivial solution $c_1=\cdots=c_k=0$. If a nontrivial combination equals zero, at least one vector can be written using the others, so the collection is dependent.

In data, dependence corresponds to exact redundancy. If one feature column is twice another, they do not supply two independent directions. A linear model cannot uniquely decide how to distribute weight between them. Near dependence, often called multicollinearity in regression, is not exact redundancy but can still make coefficients unstable.

Rank counts independent directions

The rank of a matrix is the number of linearly independent columns. It is also the number of linearly independent rows, a non-obvious but fundamental theorem. For $\mathbf{A}\in\mathbb{R}^{m\times n}$,

$$\operatorname{rank}(\mathbf{A})\le \min(m,n).$$

A matrix has full column rank if its rank equals $n$, and full row rank if its rank equals $m$. A square matrix is invertible exactly when it has full rank.

Rank provides a language for information capacity. A $100\times20$ matrix has at most 20 independent column directions. If its rank is only 7, all 20 features actually lie in a seven-dimensional linear subspace. Low-rank approximations exploit this structure to compress data and parameters. In practice, numerical rank depends on a tolerance because floating-point values that should be zero may appear as tiny nonzero singular values.

Dot products, cosine similarity, and projections

Dot products combine magnitude and direction

For vectors $\mathbf{x},\mathbf{y}\in\mathbb{R}^n$, the dot product is

$$\mathbf{x}\cdot\mathbf{y} = \mathbf{x}^T\mathbf{y} = \sum_{i=1}^n x_i y_i.$$

Geometrically,

$$\mathbf{x}\cdot\mathbf{y} = \|\mathbf{x}\|_2\|\mathbf{y}\|_2\cos\theta,$$

where $\theta$ is the angle between them. A positive dot product indicates an acute angle, zero indicates a right angle, and a negative value indicates an obtuse angle. Vectors with zero dot product are orthogonal.

Weighted sums are dot products. A linear model predicts

$$\hat y=\mathbf{w}^T\mathbf{x}+b.$$

The prediction is large when the input aligns strongly with the weight direction, adjusted by their magnitudes and the bias. This geometric view explains why a linear classifier uses a hyperplane: all points satisfying $\mathbf{w}^T\mathbf{x}+b=0$ lie on the decision boundary, and $\mathbf{w}$ is perpendicular to it.

Cosine similarity focuses on direction

Cosine similarity divides out vector lengths:

$$\operatorname{cos\_sim}(\mathbf{x},\mathbf{y}) = \frac{\mathbf{x}^T\mathbf{y}} {\|\mathbf{x}\|_2\|\mathbf{y}\|_2}.$$

For nonzero real vectors, its value lies between $-1$ and $1$. Values near $1$ mean similar direction, values near $0$ mean roughly orthogonal direction, and values near $-1$ mean opposite direction.

This makes cosine similarity useful for embeddings. A long document and a short sentence may have different embedding magnitudes but similar semantic directions. However, cosine similarity is not automatically the right metric for every embedding model. Some models encode useful information in vector magnitude, and many retrieval systems are explicitly trained for dot-product similarity. The training objective and model documentation should guide the metric choice.

Projection finds the closest component

The projection of $\mathbf{x}$ onto a nonzero vector $\mathbf{v}$ is

$$\operatorname{proj}_{\mathbf{v}}(\mathbf{x}) = \frac{\mathbf{x}^T\mathbf{v}} {\mathbf{v}^T\mathbf{v}}\mathbf{v}.$$

The coefficient tells how much of $\mathbf{x}$ points along $\mathbf{v}$. The residual

$$\mathbf{r} = \mathbf{x}-\operatorname{proj}_{\mathbf{v}}(\mathbf{x})$$

is orthogonal to $\mathbf{v}$. If $\mathbf{v}$ is already a unit vector, the formula simplifies to $(\mathbf{x}^T\mathbf{v})\mathbf{v}$.

Projection onto a subspace generalizes this idea. Least-squares regression projects the target vector $\mathbf{y}$ onto the column space of the design matrix $\mathbf{X}$. The fitted values $\hat{\mathbf{y}}=\mathbf{X}\hat{\boldsymbol{\beta}}$ are the reachable vector closest to $\mathbf{y}$. At the optimum, the residual is orthogonal to every column of $\mathbf{X}$:

$$\mathbf{X}^T(\mathbf{y}-\mathbf{X}\hat{\boldsymbol{\beta}})=\mathbf{0}.$$

This condition produces the normal equations, although stable software usually solves least squares without explicitly forming an inverse.

Eigenvalues and eigenvectors

Directions a transformation preserves

Most matrix transformations change both the length and direction of a vector. An eigenvector is a special nonzero vector whose direction is preserved:

$$\mathbf{A}\mathbf{v}=\lambda\mathbf{v}.$$

The scalar $\lambda$ is the corresponding eigenvalue. If $\lambda=3$, the transformation triples the vector along that direction. If $\lambda=-2$, it doubles the length and reverses direction. If $\lambda=0$, it collapses that direction to zero.

Eigenvectors are found by rearranging:

$$(\mathbf{A}-\lambda\mathbf{I})\mathbf{v}=\mathbf{0}.$$

For a nonzero solution to exist, $\mathbf{A}-\lambda\mathbf{I}$ must be singular, leading to the characteristic equation

$$\det(\mathbf{A}-\lambda\mathbf{I})=0.$$

The computation is useful for small examples, but intuition matters more initially. Imagine data points forming a long, tilted cloud. Some directions through the cloud show large variation; perpendicular directions show less. The covariance matrix transforms vectors in a way whose eigenvectors identify these principal directions, and whose eigenvalues measure the variance along them.

A preview of PCA

Principal component analysis starts with centered data. If $\mathbf{X}$ has observations as rows and every feature mean has been subtracted, a covariance matrix is proportional to

$$\mathbf{X}^T\mathbf{X}.$$

This matrix is symmetric and positive semidefinite, so its eigenvalues are real and nonnegative, and its eigenvectors can be chosen orthonormally. PCA orders the eigenvectors from largest eigenvalue to smallest. The first principal component is the direction of greatest variance, the second is the greatest remaining orthogonal direction, and so on.

Projecting data onto the first few eigenvectors reduces dimensionality while retaining as much variance as possible under the PCA criterion. This is useful for visualization, compression, denoising, and preprocessing. It does not guarantee that the retained variance is the most useful information for a supervised target, and it can be distorted when features use incompatible scales. Standardization may therefore be necessary before PCA.

Not every matrix has enough real eigenvectors to form a basis, and nonsymmetric matrices can have complex eigenvalues. For data analysis, singular value decomposition is often more general and numerically direct. Still, eigenvectors provide the key intuition: some directions reveal the natural action of a transformation.

A worked numeric example

Let

$$\mathbf{A}= \begin{bmatrix} 2&1\\ 1&3 \end{bmatrix}, \qquad \mathbf{x}= \begin{bmatrix} 4\\ -1 \end{bmatrix}.$$

Multiplying gives

$$\mathbf{y}=\mathbf{A}\mathbf{x} = \begin{bmatrix} 2(4)+1(-1)\\ 1(4)+3(-1) \end{bmatrix} = \begin{bmatrix} 7\\ 1 \end{bmatrix}.$$

The row interpretation says that the first output is the dot product of $\begin{bmatrix}2&1\end{bmatrix}$ with $\mathbf{x}$, while the second uses $\begin{bmatrix}1&3\end{bmatrix}$. The column interpretation says

$$\mathbf{A}\mathbf{x} = 4\begin{bmatrix}2\\1\end{bmatrix} -1\begin{bmatrix}1\\3\end{bmatrix} = \begin{bmatrix}7\\1\end{bmatrix}.$$

Now consider the squared length before and after transformation:

$$\|\mathbf{x}\|_2=\sqrt{4^2+(-1)^2}=\sqrt{17}, \qquad \|\mathbf{y}\|_2=\sqrt{7^2+1^2}=\sqrt{50}.$$

This matrix stretches $\mathbf{x}$ overall, although its effect depends on direction.

The determinant is

$$\det(\mathbf{A})=2(3)-1(1)=5,$$

so $\mathbf{A}$ is invertible. Its inverse is

$$\mathbf{A}^{-1} = \frac{1}{5} \begin{bmatrix} 3&-1\\ -1&2 \end{bmatrix}.$$

Recovering the input verifies the calculation:

$$\mathbf{A}^{-1}\mathbf{y} = \frac{1}{5} \begin{bmatrix} 3&-1\\ -1&2 \end{bmatrix} \begin{bmatrix} 7\\ 1 \end{bmatrix} = \frac{1}{5} \begin{bmatrix} 20\\ -5 \end{bmatrix} = \begin{bmatrix} 4\\ -1 \end{bmatrix}.$$

We can also compute the eigenvalues:

$$\det(\mathbf{A}-\lambda\mathbf{I}) = (2-\lambda)(3-\lambda)-1 = \lambda^2-5\lambda+5.$$

Thus

$$\lambda_{1,2}=\frac{5\pm\sqrt{5}}{2}.$$

Both are positive, which is consistent with this symmetric matrix stretching rather than reversing its eigenvector directions. These eigenvalues quantify the matrix's scaling along two perpendicular eigenvectors.

The same calculations in NumPy are concise:

import numpy as np

A = np.array([[2.0, 1.0],
              [1.0, 3.0]])
x = np.array([4.0, -1.0])

y = A @ x
x_recovered = np.linalg.solve(A, y)
eigenvalues, eigenvectors = np.linalg.eigh(A)  # for symmetric matrices

print("y:", y)
print("recovered x:", x_recovered)
print("rank:", np.linalg.matrix_rank(A))
print("eigenvalues:", eigenvalues)

The @ operator performs matrix multiplication. The elementwise expression A * x means something different: NumPy broadcasts x across the rows and multiplies corresponding entries.

Broadcasting and shapes in practice

Shapes are part of the meaning

In mathematics, a vector is often treated as a column by default. Array libraries must be more explicit. In NumPy, an array with shape (3,) is one-dimensional; it is neither a stored row matrix (1, 3) nor a column matrix (3, 1).

import numpy as np

x = np.array([1.0, 2.0, 3.0])       # shape (3,)
row = x.reshape(1, 3)                # shape (1, 3)
column = x.reshape(3, 1)             # shape (3, 1)

print(row @ column)                  # shape (1, 1)
print(column @ row)                  # shape (3, 3), an outer product
print(x @ x)                         # scalar dot product

The products differ dramatically. A row times a column produces one inner product; a column times a row produces an outer-product matrix. Relying on visual appearance alone is risky because printed one-dimensional arrays do not show orientation.

For a batch-oriented layer, choose and document a convention. A common convention is:

• X: (batch_size, input_features)
• W: (input_features, output_features)
• b: (output_features,)
• Z = X @ W + b: (batch_size, output_features)

This convention differs by a transpose from the earlier column-vector equation, but both are valid. Consistency matters more than choosing one universal orientation.

Broadcasting is implicit replication

Broadcasting allows operations on compatible shapes without explicitly copying values. If Z has shape (64, 10) and b has shape (10,), then Z + b adds the same ten-component bias to each of the 64 rows.

NumPy compares dimensions from right to left. Two dimensions are compatible when they are equal or when one of them is 1. Missing leading dimensions behave as if they were 1. Therefore:

• (64, 10) and (10,) broadcast to (64, 10).
• (64, 10) and (64, 1) broadcast to (64, 10).
• (64, 10) and (64,) are incompatible because the final dimensions 10 and 64 conflict.

Broadcasting is powerful, but a legal operation is not necessarily the intended one. Subtracting an array of shape (64, 1) from (64, 10) is valid and affects each row, whereas a mistaken reshape could silently produce a (64, 64) result. Assertions make assumptions visible:

assert X.ndim == 2
assert W.ndim == 2
assert X.shape[1] == W.shape[0]
Z = X @ W
assert Z.shape == (X.shape[0], W.shape[1])

When debugging, inspect shapes before inspecting individual numbers. Shape errors often reveal the conceptual mistake immediately.

Common mistakes and how to avoid them

Dimension mismatch

The product $(m\times n)(p\times q)$ is undefined unless $n=p$. A common error is to align dimensions that merely look related, such as batch sizes, rather than checking the output features of the first operation against the input features of the second.

Write shapes beside every variable before coding. For example,

$$\mathbf{X}_{B\times d}\mathbf{W}_{d\times h} = \mathbf{H}_{B\times h}.$$

This small habit catches many bugs before runtime.

Mixing row and column conventions

Textbooks commonly use samples as columns in some derivations and rows in data matrices elsewhere. Libraries and papers also vary. Both conventions are mathematically sound, but copying formulas across conventions without transposing weights creates confusion.

State what each axis means. Avoid saying only “a 2D tensor”; say “rows are observations and columns are features.” In code, names such as n_samples, n_features, and hidden_size communicate more than anonymous dimensions.

Confusing matrix and elementwise multiplication

Matrix multiplication combines a row with a column through a sum. Elementwise multiplication multiplies entries in matching positions. In NumPy, A @ B is matrix multiplication and A * B is elementwise multiplication. In mathematical writing, juxtaposition usually means matrix multiplication, while the Hadamard product may be written $\mathbf{A}\odot\mathbf{B}$.

The two can occasionally return arrays of the same shape, making this bug especially deceptive. Verify the intended operation, not just the output dimensions.

Assuming an inverse always exists

Only square, full-rank matrices have ordinary inverses. A feature matrix is frequently rectangular, and correlated columns can make $\mathbf{X}^T\mathbf{X}$ singular or ill-conditioned. The pseudoinverse offers a generalized solution, but it should not become a reflex that hides data problems.

Use linear solvers or least-squares routines, examine rank and conditioning, and consider regularization when appropriate. Ridge regression adds a positive term:

$$\mathbf{X}^T\mathbf{X}+\lambda\mathbf{I},$$

which improves conditioning when $\lambda>0$.

Treating floating-point equality as exact

Computer arithmetic is approximate. A theoretically orthogonal dot product may be 2e-16 rather than exactly zero. Testing value == 0 can therefore fail unexpectedly. Use tolerances such as np.isclose and np.allclose, with thresholds chosen for the scale and precision of the problem.

Similarly, numerical rank is determined relative to a tolerance. Tiny singular values should not always be interpreted as meaningful independent dimensions.

Ignoring feature scales

Euclidean distance, dot products, covariance, and PCA are sensitive to scale. If one feature is measured in thousands and another between zero and one, the first can dominate. Standardization transforms a feature using

$$z=\frac{x-\mu}{\sigma},$$

giving it mean zero and standard deviation one when $\sigma\ne0$. Standardization is not mandatory for every model, but its geometric effect must be understood. Fit scaling parameters only on training data to avoid leakage.

Losing the batch dimension

Indexing one sample may change shape from (1, d) to (d,). Later code can then broadcast or multiply differently. Use slices such as X[i:i+1] when a two-dimensional batch is required, and check framework options such as keepdims=True for reductions.

Forgetting that order matters

Because $\mathbf{A}\mathbf{B}\ne\mathbf{B}\mathbf{A}$ in general, rearranging factors is not a harmless algebraic step. The valid identity

$$(\mathbf{A}\mathbf{B})^{-1}=\mathbf{B}^{-1}\mathbf{A}^{-1}$$

also reverses order. Keep dimensions visible when manipulating expressions; invalid rearrangements frequently expose themselves through incompatible shapes.

Practice with Solver360

Reading equations builds familiarity, but calculating and visualizing transformations builds intuition. Use the free Matrix Operations Calculator to test small matrices and check hand calculations.

A productive practice sequence is:

1. Enter a $2\times2$ matrix and multiply it by several vectors. Sketch each input and output arrow.
2. Compare a diagonal scaling matrix, a rotation matrix, a shear, and a projection. Describe what each does geometrically.
3. Multiply two transformation matrices in both orders. Find a vector for which the outputs visibly differ.
4. Compute the transpose and inverse of a full-rank matrix, then verify that $\mathbf{A}^{-1}\mathbf{A}=\mathbf{I}$.
5. Try a matrix with one row equal to twice another. Observe its rank and the failure of inversion.
6. Check a manual $2\times3$ by $3\times2$ multiplication, tracking every output dot product.
7. Explore eigenvectors of a symmetric $2\times2$ matrix by confirming that their directions remain unchanged after transformation.

After using the calculator, repeat selected examples in NumPy. The calculator emphasizes structure; code emphasizes shape discipline and numerical behavior. When their answers differ, inspect orientation, operation type, and input precision before assuming either tool is wrong.

Try also constructing examples rather than only solving supplied ones. Design a matrix that reflects points across the horizontal axis. Design one that projects every vector onto the line $y=x$. Build two independent vectors and then add a third that makes the set dependent. Creating examples requires a deeper level of understanding than recognizing definitions.

Frequently asked questions

Do I need to master proofs before learning machine learning?

You do not need a full proof-based linear algebra course before training your first model. You should, however, understand definitions, dimensions, and geometric meanings well enough to predict what an operation will do. Begin with vectors, matrix multiplication, dot products, span, rank, projections, and eigenvector intuition. Implement small examples and explain the results in words.

Proofs become increasingly valuable as you study optimization, statistical learning theory, and advanced model architecture. They teach you which assumptions support a result and prevent formulas from becoming rituals. A balanced route is to learn computations and applications first, then revisit important claims—such as why least-squares residuals are orthogonal to the column space—with derivations and proofs.

What is the difference between a vector, matrix, and tensor?

A scalar has no axes, a vector has one axis, and a matrix has two axes. “Tensor” is the broader term used in machine-learning libraries for multidimensional arrays. A color image might have shape (height, width, channels), while a batch of images might have (batch, height, width, channels). Deep-learning frameworks call even scalars and matrices tensors.

In more advanced mathematics, a tensor has a transformation-based definition richer than “multidimensional array.” For everyday model implementation, shape, axis meaning, and data type are the immediate concerns. Always document what each axis represents, because two tensors with the same shape can carry entirely different semantics.

Why can a stack of linear layers be replaced by one linear layer?

Suppose two layers have no nonlinear activations:

$$\mathbf{h}=\mathbf{W}_1\mathbf{x}+\mathbf{b}_1,$$ $$\mathbf{y}=\mathbf{W}_2\mathbf{h}+\mathbf{b}_2.$$

Substitution gives

$$\mathbf{y} = \mathbf{W}_2\mathbf{W}_1\mathbf{x} +\mathbf{W}_2\mathbf{b}_1+\mathbf{b}_2.$$

Define $\mathbf{W}=\mathbf{W}_2\mathbf{W}_1$ and $\mathbf{b}=\mathbf{W}_2\mathbf{b}_1+\mathbf{b}_2$. The two layers are exactly one affine layer. Any number of affine layers can collapse in the same way. Nonlinear activations prevent this collapse and let networks bend input space into more complex decision regions. Depth without nonlinearity may still impose useful factorization constraints during optimization, but it does not expand the set of representable linear functions in the usual dense setting.

Why should I avoid explicitly computing a matrix inverse?

Computing an inverse solves a more general problem than is needed when you only want the solution to $\mathbf{A}\mathbf{x}=\mathbf{b}$. Factorization-based solvers generally require less work, introduce less numerical error, and can exploit matrix structure. For example, a symmetric positive-definite system can use Cholesky factorization.

An explicit inverse is appropriate in symbolic derivations and some small analytical discussions. In production numerical code, prefer solve, lstsq, or specialized decomposition routines. This advice does not mean inverses are unimportant; the concept explains uniqueness and reversibility. It means that the mathematically equivalent computational route can differ greatly in stability.

What does rank deficiency mean for a machine-learning dataset?

Rank deficiency means that some columns or rows are exact linear combinations of others. If a design matrix has dependent feature columns, multiple coefficient vectors may produce exactly the same predictions. For example, including temperature in both Celsius and Fahrenheit together with an intercept creates an exact linear relationship. The data cannot uniquely identify separate effects for all those columns.

Near rank deficiency is also important. Highly correlated features can lead to unstable estimates even without exact dependence. Remedies depend on the goal: remove redundant features, combine them using domain knowledge, use regularization, collect more informative data, or apply a dimensionality-reduction method. Rank is a property of the represented data, so blindly applying a numerical fix may conceal a modeling issue.

Are cosine similarity and Euclidean distance interchangeable?

Not generally. Cosine similarity considers the angle between nonzero vectors and ignores their overall magnitudes. Euclidean distance measures absolute straight-line separation and changes when either vector is rescaled. However, if every vector has been normalized to unit $L_2$ norm, the two are directly related:

$$\|\hat{\mathbf{x}}-\hat{\mathbf{y}}\|_2^2 = 2-2\hat{\mathbf{x}}^T\hat{\mathbf{y}} = 2-2\cos\theta.$$

For normalized vectors, ranking by smallest Euclidean distance is therefore equivalent to ranking by largest cosine similarity. Without normalization, this equivalence disappears. Choose a metric that matches how the embedding model was trained and what magnitude means in the application.

What is the practical difference between eigenvalue decomposition and singular value decomposition?

Eigenvalue decomposition applies directly to square matrices and expresses suitable matrices in terms of directions preserved by the transformation. Singular value decomposition applies to any $m\times n$ matrix:

$$\mathbf{A}=\mathbf{U}\boldsymbol{\Sigma}\mathbf{V}^T.$$

The right singular vectors in $\mathbf{V}$ describe input directions, the singular values describe nonnegative scaling strengths, and the left singular vectors in $\mathbf{U}$ describe corresponding output directions. SVD always exists and reveals rank, low-rank structure, and conditioning.

PCA can be computed through eigenvectors of a covariance matrix or directly using SVD of the centered data matrix. The SVD route often avoids explicitly forming $\mathbf{X}^T\mathbf{X}$ and can be numerically preferable. Eigenvectors remain useful for understanding symmetric operators, dynamics, and many theoretical results.

How much shape memorization is necessary?

Memorize principles, not every library-specific layout. For matrix multiplication, the inner dimensions must match and the outer dimensions determine the result. For a dense layer, the feature dimension of the input must meet the input dimension of the weights. For broadcasting, compare dimensions from the right and require equality or a dimension of one.

Then annotate shapes in derivations and inspect them in code. Experts do this routinely; it is not a beginner's crutch. Framework conventions can differ—for example, channel-first versus channel-last images—so deriving the expected shape from axis meaning is more reliable than memory.

Is linear algebra still important when automatic differentiation handles gradients?

Yes. Automatic differentiation computes derivatives, but it does not decide whether the model is correctly designed, whether axes align, or whether an operation is numerically sensible. Gradients are themselves vectors, matrices, and higher-dimensional arrays shaped like parameters. Backpropagation repeatedly multiplies by transposed Jacobians, even when the framework hides those matrices through efficient vector-Jacobian products.

Linear algebra helps diagnose exploding signals, redundant representations, ill-conditioned optimization, and memory costs. It also lets you read papers and understand why an implementation is equivalent to an equation. Automation reduces manual arithmetic; it does not remove the need for structural understanding.

Next steps: from transformations to PCA and gradients

The next natural topic is singular value decomposition and the full mathematics of PCA. Study centering, covariance, orthonormal bases, explained variance, and low-rank approximation. Work through a tiny centered dataset by hand: calculate its covariance matrix, find principal directions, project the observations, and reconstruct an approximation. Then repeat with SVD and connect singular values to the eigenvalues of $\mathbf{X}^T\mathbf{X}$.

Calculus provides the second major bridge. A derivative of a scalar with respect to a vector is a gradient:

$$\nabla_{\mathbf{w}}L= \begin{bmatrix} \frac{\partial L}{\partial w_1}\\ \vdots\\ \frac{\partial L}{\partial w_n} \end{bmatrix}.$$

The gradient points in the direction of steepest local increase, so gradient descent updates parameters by

$$\mathbf{w}_{t+1} = \mathbf{w}_t-\eta\nabla_{\mathbf{w}}L,$$

where $\eta$ is the learning rate. For vector-valued functions, derivatives become Jacobian matrices; for second derivatives of scalar functions, they become Hessian matrices. The chain rule then becomes a sequence of structured linear maps, which is the mathematical heart of backpropagation.

As you continue, keep connecting three views: arithmetic, geometry, and code. Arithmetic tells you how to calculate an output. Geometry explains directions, lengths, projections, and lost information. Code forces you to specify shapes, batches, precision, and stable algorithms. When all three views agree, linear algebra stops feeling like a collection of formulas and becomes a practical language for understanding machine learning.