Linear Algebra for Machine Learning

Phase 0 · Week 1–2 · ~20–25 hours

Linear algebra is the grammar of machine learning. Every dataset is a matrix. Every model parameter update is a vector operation. Every convolutional filter is a linear transformation. You have likely used these objects for a year; this file is about understanding what they are — not just syntactically, but geometrically. The goal is not to memorize formulas. The goal is to look at a matrix and see what it does to space.

Why Linear Algebra First

Before you can take a derivative of a loss function with respect to a weight matrix, you need to know what a matrix is. Before you can understand PCA, you need to understand what eigenvectors mean. Before you can read a paper that says “the attention mechanism computes scaled dot-product similarity,” you need to know what a dot product is measuring geometrically (projection + magnitude). None of this is possible without a working geometric intuition for linear algebra.

Core Concepts

1. Vectors and Geometric Intuition

A vector is a direction and a magnitude in space. In ML:

  • A training example is a vector in feature space

  • A model’s weights are a vector in parameter space

  • An embedding is a vector in a learned semantic space

The critical insight: Similarity between vectors = cosine of the angle between them. This is why dot products appear everywhere in attention mechanisms, nearest-neighbor search, and collaborative filtering.

import numpy as np

# Two document embeddings (simplified to 3D)
a = np.array([1.0, 0.5, 0.2])
b = np.array([0.9, 0.6, 0.1])

# Dot product = |a||b|cos(θ)
dot = np.dot(a, b)

# Cosine similarity (normalized dot product)
cosine_sim = dot / (np.linalg.norm(a) * np.linalg.norm(b))
print(f"Cosine similarity: {cosine_sim:.4f}")  # ~0.997 → nearly identical direction

ML connection: Transformer attention scores are Q @ K.T / sqrt(d_k). You can’t understand why that formula works until you understand what the dot product measures.

2. Matrix Multiplication as Transformation

A matrix is a function that transforms vectors. Matrix multiplication A @ B composes two transformations: first apply B, then apply A.

# Rotation matrix (rotate 45 degrees)
theta = np.pi / 4
R = np.array([[np.cos(theta), -np.sin(theta)],
              [np.sin(theta),  np.cos(theta)]])

v = np.array([1.0, 0.0])
v_rotated = R @ v
print(f"Original: {v}")
print(f"Rotated:  {v_rotated}")  # [0.707, 0.707]

ML connection: A neural network layer W @ x + b is a linear transformation followed by a translation. The entire forward pass is a chain of matrix transformations. Understanding this is why nn.Linear(in, out) is not magic — it’s a matrix W of shape (out, in).

Key properties to internalize (not memorize — understand):

  • Matrix multiplication is NOT commutative: A @ B ≠ B @ A

  • (A @ B).T = B.T @ A.T — transpose reverses order

  • Identity matrix I: A @ I = A (the “do nothing” transformation)

  • Inverse A⁻¹: A @ A⁻¹ = I (undo a transformation — not always possible)

3. Dot Products, Norms, and Projections

a = np.array([3.0, 4.0])
b = np.array([1.0, 0.0])

# L2 norm (Euclidean length)
norm_a = np.linalg.norm(a)       # 5.0
print(f"||a|| = {norm_a}")

# Projection of a onto b
proj_a_onto_b = (np.dot(a, b) / np.dot(b, b)) * b
print(f"Projection: {proj_a_onto_b}")  # [3.0, 0.0]

ML connection: Regularization penalizes the norm of weight vectors (L1: ||w||₁, L2: ||w||₂²). The choice between L1 and L2 is a choice about which geometry you want to constrain your parameters in.

4. Eigenvalues and Eigenvectors — The Geometry

An eigenvector of matrix A is a vector that doesn’t change direction when you apply A — only its magnitude changes. The eigenvalue λ tells you how much it scales.

A @ v = λ * v

# Covariance matrix of 2D data
data = np.random.randn(100, 2)
data[:, 0] *= 3  # stretch x-axis (more variance)
cov = np.cov(data.T)

eigenvalues, eigenvectors = np.linalg.eig(cov)

print(f"Eigenvalues: {eigenvalues}")
# Larger eigenvalue → direction of maximum variance
print(f"Eigenvectors:\n{eigenvectors}")
# These are the principal component directions

Why this matters for PCA:

  • The eigenvectors of the covariance matrix point in the directions of maximum variance

  • The eigenvalues tell you how much variance is in each direction

  • PCA = find the eigenvectors, project data onto them, keep only the top-k

Why this matters for neural networks:

  • The Hessian matrix (second derivatives) of a loss function has eigenvalues that tell you about the curvature of the loss landscape

  • Large eigenvalues → sharp curvature → small learning rates needed

  • This is why adaptive optimizers (Adam) work — they’re approximating information about the loss curvature

5. Singular Value Decomposition (SVD)

SVD is the most powerful decomposition in applied mathematics. Every matrix A (even non-square, even non-symmetric) decomposes as:

A = U Σ Vᵀ

Where:

  • U: left singular vectors (orthonormal, shape m×m)

  • Σ: singular values (diagonal, non-negative, sorted descending)

  • Vᵀ: right singular vectors (orthonormal, shape n×n)

# A 4x3 matrix (e.g., 4 users, 3 movie ratings)
A = np.array([[1, 0, 0],
              [0, 1, 0],
              [1, 1, 0],
              [0, 0, 1]], dtype=float)

U, sigma, Vt = np.linalg.svd(A, full_matrices=False)

print(f"U shape: {U.shape}")       # (4, 3)
print(f"Sigma: {sigma}")           # singular values
print(f"Vt shape: {Vt.shape}")     # (3, 3)

# Reconstruction (verify: A ≈ U @ np.diag(sigma) @ Vt)
A_reconstructed = U @ np.diag(sigma) @ Vt
print(f"Max reconstruction error: {np.max(np.abs(A - A_reconstructed)):.2e}")

# Truncated SVD (keep top 2 components = lossy compression)
k = 2
A_approx = U[:, :k] @ np.diag(sigma[:k]) @ Vt[:k, :]
print(f"Compression ratio: {k/3:.2f}")

ML connections for SVD:

  1. Recommendation systems: User-item matrices have low-rank structure. SVD finds latent factors.

  2. Dimensionality reduction: sklearn.decomposition.TruncatedSVD is literally this — keep top-k singular vectors.

  3. PCA: PCA on centered data IS the SVD of the data matrix (eigenvectors of covariance = right singular vectors of data).

  4. NLP: Latent Semantic Analysis is SVD on a term-document matrix.

  5. Weight initialization: Papers like “Spectral Norm Regularization” operate directly on singular values of weight matrices.

Geometric meaning: SVD says every linear transformation is a rotation + scaling + rotation. The singular values are the scaling factors. This is why SVD reveals the “essential structure” of a matrix.

6. PCA from Scratch (Core Algorithm)

This is the bridge between linear algebra theory and ML practice. Here’s the full algorithm:

def pca_from_scratch(X, n_components):
    """
    PCA implementation using eigendecomposition of covariance matrix.
    
    X: array of shape (n_samples, n_features)
    n_components: int, number of principal components to keep
    
    Returns: X_reduced (n_samples, n_components), components, explained_variance_ratio
    """
    # Step 1: Center the data (CRITICAL — never skip this)
    X_centered = X - np.mean(X, axis=0)
    
    # Step 2: Compute covariance matrix
    # Divide by (n-1) for unbiased estimate
    cov_matrix = np.cov(X_centered.T)  # shape: (n_features, n_features)
    
    # Step 3: Eigendecomposition
    eigenvalues, eigenvectors = np.linalg.eigh(cov_matrix)
    # eigh instead of eig: guaranteed real output for symmetric matrices
    
    # Step 4: Sort by eigenvalue (descending)
    sorted_idx = np.argsort(eigenvalues)[::-1]
    eigenvalues = eigenvalues[sorted_idx]
    eigenvectors = eigenvectors[:, sorted_idx]
    
    # Step 5: Select top-k components
    components = eigenvectors[:, :n_components]  # shape: (n_features, n_components)
    
    # Step 6: Project data
    X_reduced = X_centered @ components  # shape: (n_samples, n_components)
    
    # Step 7: Compute explained variance ratio
    explained_variance_ratio = eigenvalues[:n_components] / np.sum(eigenvalues)
    
    return X_reduced, components, explained_variance_ratio


# Test on iris dataset
from sklearn.datasets import load_iris
iris = load_iris()
X, y = iris.data, iris.target  # 150 samples, 4 features

X_reduced, components, evr = pca_from_scratch(X, n_components=2)
print(f"Original shape: {X.shape}")
print(f"Reduced shape: {X_reduced.shape}")
print(f"Explained variance ratio: {evr}")
print(f"Total variance captured: {evr.sum():.3f}")

# Verify against sklearn
from sklearn.decomposition import PCA
pca_sklearn = PCA(n_components=2)
X_sklearn = pca_sklearn.fit_transform(X)
# Note: sklearn may flip signs of components (both are valid)
print(f"Max abs diff from sklearn (sign-adjusted): {np.min([np.max(np.abs(X_reduced - X_sklearn)), np.max(np.abs(X_reduced + X_sklearn))]):.6f}")

What Most Learners Skip (That Kills Them Later)

They learn matrix operations but not matrix decompositions.

Knowing how to multiply matrices is table stakes. The substance is in decompositions: eigendecomposition tells you about symmetry and invariance; SVD tells you about rank and approximation; Cholesky tells you about positive definiteness (which appears in covariance matrices and kernel methods). If you’ve never heard of Cholesky decomposition, you will be confused when you encounter Gaussian Processes.

Second blind spot: the difference between eig and eigh. For symmetric matrices (covariance matrices ARE symmetric), np.linalg.eigh is both numerically more stable and guaranteed to return real eigenvalues. Using np.linalg.eig on a covariance matrix is a subtle bug that will produce complex-valued components and confuse you.

Third blind spot: not centering data before PCA. This is the single most common implementation error. Uncentered PCA finds components that explain variance including the mean offset — which means your first principal component often just points toward the data mean. Always center first.

Practice Problems with Acceptance Criteria

Problem 1 — Matrix Properties Implement a function that checks whether a matrix is: (a) symmetric, (b) orthogonal, (c) positive semi-definite. Test it on: a covariance matrix, a rotation matrix, and a random matrix.

  • Acceptance criteria: All three checks work correctly on all three test cases. Explain in a comment why PSD matters for covariance matrices.

Problem 2 — Eigendecomposition Verification For a 3×3 symmetric matrix you create, verify the eigendecomposition equation A @ V = V @ D where V contains eigenvectors and D is diagonal with eigenvalues. Verify with np.allclose.

  • Acceptance criteria: np.allclose(A @ V, V @ D, atol=1e-10) returns True.

Problem 3 — SVD Image Compression Load any grayscale image as a numpy array. Apply truncated SVD with k=5, 20, 50, 100 components. Plot the reconstruction quality vs compression ratio. What k gives you ~95% of variance explained?

  • Acceptance criteria: Working plot with labeled axes, correct variance calculation, stated optimal k value.

Problem 4 — PCA from Scratch (Phase Project Warm-up) Implement pca_from_scratch (code above), then verify it matches sklearn.decomposition.PCA on the iris dataset within floating point tolerance.

  • Acceptance criteria: np.allclose passes (after sign adjustment). Explained variance ratios match to 4 decimal places.

Resource Schedule

Resource

What to Cover

Time Estimate

3B1B “Essence of Linear Algebra”

All 16 chapters (21 hrs total, 13 hrs video)

8–10 hrs (watch + pause + think)

Gilbert Strang MIT 18.06 OCW

Lectures 1–14 (vectors → eigenvalues)

10–12 hrs

Gilbert Strang MIT 18.06 OCW

Lectures 15–22 (eigenvalues → SVD)

8–10 hrs

Numpy practice (problems above)

Implement all 4 problems

4–5 hrs

Total

~30–37 hrs

Why 3B1B before Strang: 3B1B gives you the geometric intuition in ~10 hours. Strang then feels like filling in the rigorous details of something you already half-know rather than encountering it cold. The reverse order (Strang first) works but is harder.

Quick Reference: Key Formulas

Concept

Formula

ML Where It Appears

Dot product

a·b = Σᵢ aᵢbᵢ = |a||b|cos(θ)

Attention, similarity, projection

L2 norm

|v|₂ = √(Σᵢ vᵢ²)

L2 regularization, distance metrics

Matrix multiply

(AB)ᵢⱼ = Σₖ Aᵢₖ Bₖⱼ

Forward pass, linear layers

Eigendecomposition

Av = λv

PCA, spectral methods

SVD

A = UΣVᵀ

Dimensionality reduction, rec systems

Covariance matrix

Σ = (1/n) Xᵀ X (centered)

PCA, Gaussian distributions

Return to README.md · Next: 02_calculus_and_optimization.md