🎯 What You'll Learn
- Represent scalars, vectors, and matrices as NumPy arrays and understand their ML meanings
- Perform vector operations: addition, scalar multiplication, and dot products using
np.dot()and@ - Multiply matrices correctly — including shape rules and why order matters — with
@ - Compute vector norms (L1 and L2) and understand their connection to regularization
- Understand the matrix inverse, identity matrix, and eigenvectors conceptually and in NumPy
This lesson closes Phase 1 by teaching the math that every later phase runs on. Along the way, it names things you haven't studied yet — linear regression (Lesson 11), neural networks (Phase 4), PCA (Phase 3). That's deliberate: treat each mention as a signpost, not something you're expected to already know. Every one of those systems, when you reach it, will turn out to be a handful of the operations covered right here.
1 Why Linear Algebra for ML?
Every dataset you work with is a matrix. If you have 1000 customers with 15 features each, you have a 1000×15 matrix X. And remember the weighted-sum grading example from Lesson 02? ML models make predictions the same way: a model is essentially a learned vector of 15 weights, β, and predicting for a new customer is one dot product: ŷ = x · β. (You'll build your first such model in Lesson 11.)
Here's a concrete mapping of ML concepts to linear algebra objects — the right column is a preview of the whole course, so don't worry if most entries are unfamiliar:
| ML Concept | Linear Algebra Object | NumPy Shape |
|---|---|---|
| Dataset (n samples, p features) | Matrix X | (n, p) |
| Model weights / coefficients | Vector β | (p,) or (p, 1) |
| Target labels | Vector y | (n,) |
| Single prediction | Dot product x · β | scalar |
| Neural network layer | Matrix multiplication + bias | Z = X @ W + b |
| Word embeddings (NLP) | Vectors in high-dim space | (vocab_size, embed_dim) |
| PCA — principal components | Eigenvectors of covariance matrix | (p, p) |
You don't need a full linear algebra course before doing ML — but understanding the core operations (dot product, matrix multiplication, transpose) and their meanings will make every algorithm make intuitive sense rather than being a black box.
2 Scalars, Vectors & Matrices
The building blocks of linear algebra correspond directly to NumPy array dimensions:
Scalars (0D)
A single number. In ML: a single prediction, a loss value, a learning rate. Shape: ().
import numpy as np
learning_rate = np.float32(0.001)
loss_value = np.float32(0.342)
print(f"Scalar: {learning_rate}, shape: {np.array(learning_rate).shape}")Vectors (1D)
An ordered list of numbers. In ML: features of a single data point, model weights, word embeddings. Shape: (n,).
import numpy as np
# Feature vector: age, income, years_employed, credit_score
customer = np.array([34, 75000, 7, 720])
print(f"Feature vector: {customer}")
print(f"Shape: {customer.shape}") # (4,)
print(f"Dimensions: {customer.ndim}") # 1
# Weight vector for a linear model
weights = np.array([0.02, 0.001, 0.15, 0.003])
print(f"\nWeight vector: {weights}")
print(f"Shape: {weights.shape}") # (4,)Matrices (2D)
A 2D grid of numbers arranged in rows and columns. In ML: entire datasets, neural network weight matrices. Shape: (rows, cols).
import numpy as np
# Dataset matrix: 4 samples, 3 features
X = np.array([
[1.2, 3.4, 0.5], # sample 0
[2.1, 1.0, 4.3], # sample 1
[0.8, 2.9, 1.1], # sample 2
[3.0, 0.5, 2.7] # sample 3
])
print(f"Dataset matrix X:")
print(X)
print(f"Shape: {X.shape}") # (4, 3) — 4 rows (samples), 3 columns (features)
print(f"n_samples: {X.shape[0]}")
print(f"n_features: {X.shape[1]}")
# Neural network weight matrix: 3 inputs → 2 neurons
W = np.random.randn(3, 2)
print(f"\nWeight matrix W shape: {W.shape}") # (3, 2)3D and Higher (Tensors)
3D and 4D arrays are called tensors. In ML: sequences (batch × timesteps × features), images (batch × height × width × channels).
import numpy as np
# Batch of 32 sequences, each 10 timesteps, 64 features (e.g., for an LSTM)
sequences = np.zeros((32, 10, 64))
print(f"Sequences shape: {sequences.shape}") # (32, 10, 64)
# Batch of 16 color images, 224x224 pixels, 3 channels
images = np.zeros((16, 224, 224, 3))
print(f"Image batch shape: {images.shape}") # (16, 224, 224, 3)3 Vector Operations
Vectors support element-wise operations and the dot product. Understanding the dot product geometrically is key to understanding similarity measures and neural networks.
Element-Wise Operations
import numpy as np
a = np.array([1, 2, 3, 4])
b = np.array([5, 6, 7, 8])
print("Element-wise addition: ", a + b) # [6, 8, 10, 12]
print("Element-wise subtraction: ", a - b) # [-4, -4, -4, -4]
print("Element-wise multiplication: ", a * b) # [5, 12, 21, 32]
print("Scalar multiplication: ", 3 * a) # [3, 6, 9, 12]Geometrically, vector addition is "tip-to-tail": place the second vector's tail at the first vector's tip, and the resultant vector runs from the origin to the new tip.
Tip-to-tail vector addition: slide vector b so its tail sits at the tip of vector a. The resultant a+b is the vector from the origin directly to the tip of b — equivalent to adding the vectors' components element-wise.
Dot Product
The dot product of two vectors a and b is defined as:
a · b = Σ aᵢ × bᵢ = |a| × |b| × cos(θ)
Two key properties emerge from the geometric form:
- If vectors point in the same direction (θ=0°): dot product is maximized (= |a||b|)
- If vectors are perpendicular (θ=90°): dot product = 0
- If vectors point in opposite directions (θ=180°): dot product is negative (= −|a||b|)
This is why cosine similarity (the normalized dot product) measures semantic similarity between word embeddings — words used in similar contexts have embeddings that point in similar directions.
Drag the slider to rotate vector b around the origin and watch the dot product and cosine similarity update live. Vector a stays fixed along the x-axis with |a|=4; vector b has |b|=3. Notice how both quantities are maximal when the vectors are aligned (θ=0°), cross zero exactly at θ=90° (perpendicular), and go negative once θ exceeds 90°:
θ = 45° — vectors point in broadly similar directions; dot product and cosine similarity are both positive.
import numpy as np
a = np.array([1, 2, 3])
b = np.array([4, 5, 6])
# Three equivalent ways to compute the dot product
dot1 = np.dot(a, b) # NumPy function
dot2 = a @ b # @ operator (preferred in modern code)
dot3 = np.sum(a * b) # manual: sum of element-wise products
print(f"np.dot(a, b): {dot1}") # 1*4 + 2*5 + 3*6 = 4 + 10 + 18 = 32
print(f"a @ b: {dot2}") # 32
print(f"manual: {dot3}") # 32
# ── ML use case: linear model prediction ──
# A weighted sum, exactly like the course-grade example in Lesson 02 —
# except here a model has LEARNED the weights instead of being given them
x = np.array([0.5, 1.2, -0.3]) # one customer's features (standardized)
weights = np.array([0.8, 1.5, 0.4]) # weights the model learned from data
bias = 0.1 # plus a fixed starting offset ("bias")
prediction = x @ weights + bias # the fundamental linear model operation
print(f"\nLinear model prediction: {prediction:.4f}")
# ── Cosine similarity: measure similarity between vectors ──
def cosine_similarity(v1, v2):
return np.dot(v1, v2) / (np.linalg.norm(v1) * np.linalg.norm(v2))
# Suppose each word is represented as a vector of numbers, where
# similar words get similar vectors (in Phase 5 you'll learn how
# such "word embeddings" are created — here, take them as given)
king = np.array([0.5, 0.7, 0.2])
queen = np.array([0.4, 0.8, 0.3])
apple = np.array([-0.6, 0.1, 0.9])
print(f"\ncos(king, queen): {cosine_similarity(king, queen):.4f}") # high: ~0.99
print(f"cos(king, apple): {cosine_similarity(king, apple):.4f}") # low: ~0.154 Matrix Multiplication
Matrix multiplication is the operation at the heart of every neural network, every linear model prediction on a batch of samples, and every PCA computation. It is not element-wise — it combines rows of the first matrix with columns of the second.
The Shape Rule
A matrix of shape (m × n) can be multiplied by a matrix of shape (n × p). The result has shape (m × p). The inner dimensions must match. If they don't, you get a shape error.
Each output cell is the dot product of one row from the first matrix and one column from the second. Using the example below (A is 2×3, B is 3×2): the highlighted row of A and highlighted column of B combine to produce a single output cell, C[0,0] = 1·7 + 2·9 + 3·11 = 58.
Each output cell is a dot product: row 0 of A (green) combines with column 0 of B (green) to produce the single output cell C[0,0] (amber). Sliding this row/column pair across A and B fills in every cell of C — row i of A always pairs with column j of B to produce C[i, j].
import numpy as np
# ── Simple matrix multiplication ──
A = np.array([[1, 2, 3],
[4, 5, 6]]) # shape (2, 3)
B = np.array([[7, 8],
[9, 10],
[11, 12]]) # shape (3, 2)
# A is (2×3), B is (3×2) → inner dims match (3=3) → result is (2×2)
C = A @ B
print("A @ B =")
print(C)
# [[ 58 64]
# [139 154]]
# Verify: C[0,0] = 1*7 + 2*9 + 3*11 = 7 + 18 + 33 = 58 ✓
print(f"\nC shape: {C.shape}") # (2, 2)
# ── ML use case: batch prediction ──
# X: 5 samples × 4 features
X = np.random.randn(5, 4)
# W: one weight column per output — this maps 4 features to 3 outputs
W = np.random.randn(4, 3)
# b: 3 fixed offsets, one per output
b = np.zeros(3)
# All 5 predictions at once — one row per sample, one column per output
Z = X @ W + b
print(f"\nX shape: {X.shape}") # (5, 4)
print(f"W shape: {W.shape}") # (4, 3)
print(f"Z shape: {Z.shape}") # (5, 3) — 5 samples, 3 output values each
# ── NOT commutative ──
A2 = np.random.randn(3, 2)
B2 = np.random.randn(2, 3)
AB = A2 @ B2 # (3×2) @ (2×3) = (3×3)
BA = B2 @ A2 # (2×3) @ (3×2) = (2×2) — completely different shape!
print(f"\nA @ B shape: {AB.shape}") # (3, 3)
print(f"B @ A shape: {BA.shape}") # (2, 2)
print("Matrix multiplication is NOT commutative: AB ≠ BA in general")When you get ValueError: matmul: Input operand 1 has a mismatch in its core dimension 0, it means the inner dimensions of your matrices don't match. Always print(X.shape, W.shape) before a matrix multiplication. The rule: (m, n) @ (n, p) → (m, p). The bold numbers must be equal.
5 Transpose
The transpose of a matrix flips its rows and columns. A matrix of shape (m, n) becomes shape (n, m) after transposing. In NumPy: A.T or A.transpose().
import numpy as np
A = np.array([[1, 2, 3],
[4, 5, 6]]) # shape (2, 3)
A_T = A.T
print("A:")
print(A)
# [[1 2 3]
# [4 5 6]]
print("\nA.T:")
print(A_T)
# [[1 4]
# [2 5]
# [3 6]]
print(f"A shape: {A.shape}, A.T shape: {A_T.shape}") # (2,3), (3,2)
# ── Why transpose matters in ML ──
# 1. Computing the covariance matrix: (1/n) * X.T @ X
# X shape: (n_samples, n_features)
X = np.random.randn(100, 5) # 100 samples, 5 features
cov_matrix = (1 / (len(X) - 1)) * X.T @ X
print(f"\nCovariance matrix shape: {cov_matrix.shape}") # (5, 5)
# This (5x5) matrix encodes how each feature varies with every other feature
# 2. The normal equations in linear regression: beta = inv(X.T @ X) @ X.T @ y
# (covered in Section 7)
# 3. Dot product via transpose: a.T @ b = a · b (for column vectors)
a = np.array([[1], [2], [3]]) # column vector, shape (3,1)
b = np.array([[4], [5], [6]]) # column vector, shape (3,1)
print(f"\na.T @ b = {(a.T @ b)[0,0]}") # 32 — same as np.dot(a.flatten(), b.flatten())6 Norms: Measuring Vector Length
A norm is a function that assigns a length (magnitude) to a vector. Different norms measure "size" differently, and each has specific ML applications.
L2 Norm (Euclidean Norm)
The L2 norm is the familiar straight-line distance from the origin: ||v||₂ = √(v₁² + v₂² + ... + vₙ²). This is the default when people say "the length of a vector." Where you'll use it: in Lesson 21, penalizing a model whose weight vector has a large L2 norm is how you keep it from over-relying on any feature (a technique called Ridge regularization).
L1 Norm (Manhattan Norm)
The L1 norm sums the absolute values: ||v||₁ = |v₁| + |v₂| + ... + |vₙ| — the distance a taxi drives on a street grid rather than as the crow flies. Where you'll use it: penalizing a large L1 norm (Lesson 21's Lasso regularization) has the interesting side effect of pushing some weights to exactly zero.
import numpy as np
v = np.array([3.0, -4.0, 0.0, 2.0])
# L2 norm: Euclidean distance from origin
l2 = np.linalg.norm(v) # default is L2
l2_manual = np.sqrt(np.sum(v**2))
print(f"L2 norm: {l2:.4f} (manual: {l2_manual:.4f})") # sqrt(9+16+0+4) = sqrt(29) ≈ 5.385
# L1 norm: sum of absolute values
l1 = np.linalg.norm(v, ord=1)
l1_manual = np.sum(np.abs(v))
print(f"L1 norm: {l1:.4f} (manual: {l1_manual:.4f})") # 3+4+0+2 = 9.0
# L∞ norm: maximum absolute value
linf = np.linalg.norm(v, ord=np.inf)
print(f"L∞ norm: {linf:.4f}") # max(3,4,0,2) = 4.0
# ── Unit vector: normalize to length 1 ──
v_unit = v / np.linalg.norm(v)
print(f"\nUnit vector: {v_unit.round(4)}")
print(f"L2 norm of unit vector: {np.linalg.norm(v_unit):.6f}") # 1.000000
# ── Preview of Lesson 21: norms as "weight penalties" ──
# Imagine a model learned these 5 weights; the norms measure
# how "big" the weight vector is overall
beta = np.array([1.5, 0.001, 2.3, 0.002, 1.1])
l1_penalty = 0.1 * np.linalg.norm(beta, ord=1) # L1-based penalty
l2_penalty = 0.1 * np.linalg.norm(beta)**2 # L2-based penalty (squared)
print(f"\nModel weights: {beta}")
print(f"L1 penalty (λ=0.1): {l1_penalty:.4f}")
print(f"L2 penalty (λ=0.1): {l2_penalty:.4f}")
# In Lesson 21 you'll add penalties like these to a model's training
# objective to keep its weights small — nothing more than these normsFile this away for later: penalizing the L1 norm tends to push many weights to exactly 0 (the model effectively ignores those features), while penalizing the L2 norm shrinks all weights proportionally without zeroing them. When Lesson 21 introduces Lasso, Ridge, and ElasticNet, this one geometric fact is the entire intuition behind them.
7 Identity Matrix & Matrix Inverse
The identity matrix I is the matrix equivalent of the number 1: multiplying any matrix A by I gives back A. It's a square matrix with 1s on the main diagonal and 0s everywhere else.
The matrix inverse A⁻¹ satisfies A × A⁻¹ = I. It's the matrix equivalent of dividing by A. Only square matrices can have inverses, and not all square matrices are invertible (singular matrices have no inverse).
import numpy as np
# ── Identity matrix ──
I3 = np.eye(3)
print("3×3 Identity matrix:")
print(I3)
# [[1. 0. 0.]
# [0. 1. 0.]
# [0. 0. 1.]]
A = np.array([[3.0, 1.0],
[2.0, 4.0]])
print("\nA @ I = A:")
print(A @ np.eye(2)) # same as A
# ── Matrix inverse ──
A_inv = np.linalg.inv(A)
print("\nA inverse:")
print(A_inv.round(4))
# Verify: A @ A_inv ≈ Identity
print("\nA @ A_inv ≈ I:")
print((A @ A_inv).round(10))
# ── Solving a linear system Ax = b ──
# This is more numerically stable than computing inv(A) explicitly
A_sys = np.array([[2.0, 1.0, -1.0],
[-3.0, -1.0, 2.0],
[-2.0, 1.0, 2.0]])
b_rhs = np.array([8.0, -11.0, -3.0])
x_solution = np.linalg.solve(A_sys, b_rhs) # preferred over inv(A) @ b
print(f"\nSolution to Ax = b: x = {x_solution}")
print(f"Verification A @ x: {(A_sys @ x_solution).round(6)}") # should equal b
# ── Singular matrix (no inverse) ──
singular = np.array([[1.0, 2.0],
[2.0, 4.0]]) # row 2 = 2 * row 1 (linearly dependent)
try:
np.linalg.inv(singular)
except np.linalg.LinAlgError as e:
print(f"\nSingular matrix error: {e}")
# Use np.linalg.pinv() for a pseudo-inverse when the matrix may be singularComputing np.linalg.inv(A) @ b is numerically less stable than np.linalg.solve(A, b). The solve function uses LU decomposition internally, which avoids the numerical errors that accumulate when explicitly inverting large or ill-conditioned matrices. In production ML code, always prefer solve when you need to solve Ax=b.
8 Eigenvectors & Eigenvalues (Conceptual)
An eigenvector of a matrix A is a special non-zero vector v that, when multiplied by A, only gets scaled — it doesn't change direction. The scaling factor is called the eigenvalue λ. Formally:
Av = λv
Most vectors get both rotated and scaled when multiplied by a matrix. Eigenvectors are the special "principal axes" that only get scaled. This geometric specialness is exactly why they're useful — here are three famous places they show up (all previews; none required yet):
- PCA (Phase 3, Lesson 31): A technique for compressing many correlated features into a few informative ones. Its "principal components" are literally eigenvectors, and the eigenvalues say how much information each one carries.
- PageRank (Google's original algorithm): Finding the most important web pages is equivalent to finding the dominant eigenvector of the web's link matrix.
- Advanced deep learning: Some neural networks that operate on graphs (networks of connected nodes) are built on eigenvectors of the graph's connection matrix.
import numpy as np
# ── Compute eigenvalues and eigenvectors ──
A = np.array([[4.0, 2.0],
[1.0, 3.0]])
eigenvalues, eigenvectors = np.linalg.eig(A)
print("Eigenvalues: ", eigenvalues) # [5. 2.]
print("Eigenvectors (columns):")
print(eigenvectors)
# [[ 0.894 -0.707]
# [ 0.447 0.707]]
# Verify: A @ v = λ * v for the first eigenvector
v1 = eigenvectors[:, 0] # first eigenvector (column)
lam1 = eigenvalues[0] # first eigenvalue
print(f"\nA @ v1: {A @ v1}")
print(f"λ1 * v1: {lam1 * v1}")
print(f"Are they equal? {np.allclose(A @ v1, lam1 * v1)}") # True
# ── PCA connection: eigenvectors of the covariance matrix ──
np.random.seed(42)
X = np.random.randn(200, 3)
X[:, 1] = X[:, 0] * 0.8 + np.random.randn(200) * 0.3 # feature 1 correlates with 0
# Covariance matrix (symmetric)
X_centered = X - X.mean(axis=0)
cov = (X_centered.T @ X_centered) / (len(X) - 1)
eigenvalues_cov, eigenvectors_cov = np.linalg.eig(cov)
# Sort by descending eigenvalue (most important component first)
idx = np.argsort(eigenvalues_cov)[::-1]
eigenvalues_sorted = eigenvalues_cov[idx]
eigenvectors_sorted = eigenvectors_cov[:, idx]
print(f"\nVariance explained by each component:")
for i, ev in enumerate(eigenvalues_sorted):
pct = ev / eigenvalues_sorted.sum() * 100
print(f" PC{i+1}: {pct:.1f}%")
# Project data onto first 2 principal components (dimensionality reduction)
X_pca = X_centered @ eigenvectors_sorted[:, :2]
print(f"\nOriginal shape: {X.shape} → After PCA (2 components): {X_pca.shape}")When you study PCA properly in Phase 3, you'll use a ready-made tool (sklearn.decomposition.PCA) rather than computing eigenvectors by hand. The example above is the math running inside that tool. Keep one sentence from it: the first principal component is the eigenvector with the largest eigenvalue — the direction along which the data varies the most.
Real-World Spotlight: Linear Algebra in Production ML Systems
Here's how the operations from this lesson map directly to three real production ML components: linear regression (Lesson 11), a neural network forward pass (Phase 4), and word embedding similarity (Phase 5). You haven't studied any of these systems yet — that's the point. Read the code line by line and notice that every single operation is one you now know. When you meet these systems for real, the math will already be an old friend.
import numpy as np
np.random.seed(42)
# ════════════════════════════════════════════════
# 1. LINEAR REGRESSION — Ordinary Least Squares
# ════════════════════════════════════════════════
# The closed-form OLS solution: β = (XᵀX)⁻¹ Xᵀy
# This is matrix math — no loops, computes all weights at once
n_samples, n_features = 200, 4
# Simulate: house price depends on area, rooms, age, distance_to_city
X_raw = np.random.randn(n_samples, n_features)
true_weights = np.array([3.5, 1.2, -0.8, -1.5]) # true relationship
y = X_raw @ true_weights + np.random.randn(n_samples) * 0.5
# Add bias column (column of ones)
X = np.column_stack([np.ones(n_samples), X_raw]) # shape (200, 5)
# OLS solution: β = (XᵀX)⁻¹ Xᵀy
XtX = X.T @ X # (5, 5) — called the Gram matrix
Xty = X.T @ y # (5,)
beta = np.linalg.solve(XtX, Xty) # more stable than inv(XtX) @ Xty
print("=== OLS Linear Regression ===")
print(f"True weights: {true_weights}")
print(f"Recovered: {beta[1:].round(3)}") # skip bias term
print(f"Bias term: {beta[0]:.3f}")
# Predict on training data
y_pred = X @ beta
ss_res = np.sum((y - y_pred) ** 2)
ss_tot = np.sum((y - y.mean()) ** 2)
r2 = 1 - ss_res / ss_tot
print(f"R² score: {r2:.4f}")
# ════════════════════════════════════════════════
# 2. NEURAL NETWORK FORWARD PASS (2-layer network)
# ════════════════════════════════════════════════
# A neural network (Phase 4) is layers of exactly the batch-prediction
# pattern from Section 4: matrix multiply, add offsets, then apply a
# simple element-wise function between layers. That's all this is.
def relu(x):
return np.maximum(0, x) # zero out negatives (Lesson 02's np.where trick)
def sigmoid(x):
return 1 / (1 + np.exp(-x)) # squash any number into (0, 1)
# Dimensions: 4 inputs → 8 hidden neurons → 1 output
W1 = np.random.randn(4, 8) * 0.1 # shape (4, 8)
b1 = np.zeros(8) # shape (8,)
W2 = np.random.randn(8, 1) * 0.1 # shape (8, 1)
b2 = np.zeros(1) # shape (1,)
# Forward pass for a batch of 32 samples
batch = np.random.randn(32, 4) # 32 samples, 4 features
Z1 = batch @ W1 + b1 # (32, 4) @ (4, 8) = (32, 8)
A1 = relu(Z1) # activation — shape still (32, 8)
Z2 = A1 @ W2 + b2 # (32, 8) @ (8, 1) = (32, 1)
A2 = sigmoid(Z2) # output probabilities — shape (32, 1)
print("\n=== Neural Network Forward Pass ===")
print(f"Input batch: {batch.shape}")
print(f"After layer 1: {A1.shape} (Z1={Z1.shape}, ReLU applied)")
print(f"Output: {A2.shape} (sigmoid probabilities)")
print(f"Predictions: min={A2.min():.3f}, max={A2.max():.3f}")
# ════════════════════════════════════════════════
# 3. WORD EMBEDDINGS — Dot Product as Similarity
# ════════════════════════════════════════════════
# Simulated word embedding matrix (50-dim embeddings, 6 words)
vocab = ['king', 'queen', 'man', 'woman', 'prince', 'automobile']
embeddings = np.random.randn(6, 50)
# Make 'king' and 'queen' semantically similar (share direction)
embeddings[1] = embeddings[0] + np.random.randn(50) * 0.1
# Make 'prince' partially similar to 'king'
embeddings[4] = embeddings[0] * 0.7 + np.random.randn(50) * 0.3
def cosine_sim(a, b):
return np.dot(a, b) / (np.linalg.norm(a) * np.linalg.norm(b))
print("\n=== Word Embedding Similarities ===")
query = embeddings[0] # 'king'
for i, word in enumerate(vocab):
sim = cosine_sim(query, embeddings[i])
print(f" cos(king, {word:12s}) = {sim:.4f}")These three examples cover the most important linear algebra operations you'll encounter in production: the OLS solution for linear regression uses X.T @ X and np.linalg.solve(); the neural network forward pass is a sequence of @ W + b operations; and word embedding similarity is just a normalized dot product. Every algorithm in Phase 2 and beyond builds on exactly these operations.
✍️ Practice Exercise
Implement the following linear algebra operations from scratch using NumPy (no scikit-learn for these exercises):
- Create two 3D vectors u = [1, 2, 3] and v = [4, 5, 6]. Compute: their dot product, the L2 norm of each, and their cosine similarity. Interpret the cosine similarity value.
- Create a (4×3) matrix A and a (3×5) matrix B using
np.random.randn. Compute their product C = A @ B. Verify the shape. What happens if you try B @ A? - Solve the following linear system using
np.linalg.solve():2x + y = 5 x + 3y = 10
Verify your solution by plugging back in. - Generate a random 100×5 dataset X. Center it (subtract the mean of each column). Compute the 5×5 covariance matrix as
X_c.T @ X_c / (n-1). Compute its eigenvalues. What percentage of total variance does the first principal component explain?
▶ Show Solution
import numpy as np
# Task 1: vector operations
u = np.array([1, 2, 3], dtype=float)
v = np.array([4, 5, 6], dtype=float)
dot = u @ v
norm_u = np.linalg.norm(u)
norm_v = np.linalg.norm(v)
cos_sim = dot / (norm_u * norm_v)
print(f"Dot product u·v: {dot:.2f}")
print(f"||u||₂: {norm_u:.4f}")
print(f"||v||₂: {norm_v:.4f}")
print(f"Cosine similarity: {cos_sim:.4f}") # ~0.97 — very similar direction
# Task 2: matrix multiplication
np.random.seed(1)
A = np.random.randn(4, 3)
B = np.random.randn(3, 5)
C = A @ B
print(f"\nA@B shape: {C.shape}") # (4, 5)
try:
B @ A # (3,5) @ (4,3) — inner dims 5≠4, fails
except ValueError as e:
print(f"B@A error: {e}")
# Task 3: solve linear system
# 2x + y = 5
# x + 3y = 10
A_sys = np.array([[2, 1], [1, 3]], dtype=float)
b_rhs = np.array([5, 10], dtype=float)
sol = np.linalg.solve(A_sys, b_rhs)
print(f"\nSolution: x={sol[0]:.4f}, y={sol[1]:.4f}") # x=1, y=3
print(f"Verification A@x: {A_sys @ sol}") # should be [5, 10]
# Task 4: PCA via covariance matrix
np.random.seed(0)
X = np.random.randn(100, 5)
X[:, 1] = X[:, 0] * 0.9 + np.random.randn(100) * 0.2
X_c = X - X.mean(axis=0) # center
cov = X_c.T @ X_c / (len(X) - 1) # (5, 5) covariance
eigenvalues, _ = np.linalg.eig(cov)
eigenvalues = np.sort(eigenvalues.real)[::-1] # sort descending
total_var = eigenvalues.sum()
for i, ev in enumerate(eigenvalues):
print(f"PC{i+1}: {ev/total_var:.1%} variance explained")📚 Primary Source for This Lesson
NumPy Linear Algebra — official NumPy documentation
Complete reference for all np.linalg functions: matrix operations, decompositions, norms, and solvers. Equally essential: Chapter 2 of the Deep Learning Book (Goodfellow et al.) — the definitive treatment of linear algebra for ML, freely available online.