Unit 1.8 · about 6 h
58 / 324Batched shapes, operation costs and conditioning
Shapes determine what is computed and how much work it takes. Conditioning determines how much information small input errors can destroy even when the algorithm is stable.
- derive
- measure
You can now
- Trace batched products, transposes, reductions and broadcasting with explicit axis names.
- Derive parameter, activation and arithmetic counts for dense maps and distinguish factorization from solve costs.
- Prove the fixed-\(A\) perturbation bound and condition-number squaring from the SVD.
- Measure relative residual and solution error on a controlled failure case without confusing conditioning with algorithmic stability.
Covers
- Batched shapes, broadcasting, parameter counts and operation costs
- Condition number, residual error, solution error and numerical stability
Study, in this order
6 steps
- Broadcast-compatible axesNumPy broadcasting and matmul shape rules
General broadcasting rules and Broadcastable arrays
Depth · Stop before the vector-quantization application.
- Batched matrix productsNumPy broadcasting and matmul shape rules
Notes: two-dimensional multiplication, stacked matrices, one-dimensional promotion, and the comparison with
dot - Count products and QR workIntroduction to Applied Linear Algebra
Chapter 10, section 10.1, Complexity of matrix multiplication and Complexity of matrix triple product, printed pp. 182 to 183; Chapter 12, section 12.3, Complexity, p. 232
Depth · Use dense arithmetic counts. Leave specialized sparse products and do not interpret these counts as measured runtimes.
- Factor once, solve repeatedlyFundamentals of Numerical Computation: conditioning and operation counts
Python edition, section 2.5.3, Solution of linear systems: Lemma 2.5.1, equations (2.5.7) to (2.5.8), and Theorem 2.5.1
Depth · Derive the loop sums for triangular solves and LU. Skip the unseeded timing demonstrations and leave general asymptotic-definition proofs to Algorithms.
- Condition number and the residual-error gapFundamentals of Numerical Computation: conditioning and operation counts
Python edition, section 2.8, equations (2.8.1) to (2.8.6), Definition 2.8.1, then section 2.8.2, equations (2.8.9) to (2.8.11)
Depth · Prove the finite bound for perturbing b with A fixed. Skip the unseeded Hilbert example and the first-order matrix-perturbation formula (2.8.7); a finite bound for simultaneous perturbations is outside this unit.
- Norms and solver difficultyMIT 18.065 - Applied matrix bridge
Lecture 8, Norms of Vectors and Matrices; Lecture 10, Survey of Difficulties with \(Ax=b\)
Depth · Focus on the induced 2-norm, conditioning and loss of accuracy. Leave iterative and large-data solution methods.
Focus on
- Write the contracted axes separately from broadcast batch axes. For a dense \(d\)-to-\(h\)-to-\(p\) map, shared parameters with biases total \(dh+h+hp+p\); retaining input, hidden and output for \(N\) vectors stores \(N(d+h+p)\) scalar activations.
- The standard dot-product algorithm costs \(mn(2k-1)\) operations for an \(m\) by \(k\) times \(k\) by \(n\) product. Modified Gram-Schmidt QR costs about \(2mn^2\) for \(m\ge n\); optimized LU costs about \((2/3)n^3\) and each triangular solve about \(n^2\). These are algorithm-specific counts, not timing claims.
- For invertible square \(A\), \(\kappa_2(A)=\sigma_{\max}/\sigma_{\min}\). From \(\Delta x=A^{-1}\Delta b\) and \(\|b\|\le\|A\|\|x\|\) derive \(\|\Delta x\|/\|x\|\le\kappa(A)\|\Delta b\|/\|b\|\) for nonzero \(b\), using compatible induced norms.
- For full-column-rank rectangular \(A\), define \(\kappa_2(A)=\sigma_{\max}/\sigma_{\min}\) over its \(n\) positive singular values. The SVD gives eigenvalues \(\sigma_i^2\) for \(A^\top A\), hence \(\kappa_2(A^\top A)=\kappa_2(A)^2\). Rank deficiency makes the Gram matrix singular.
Reproduce closed-book
float64 X=np.arange(24,dtype=np.float64).reshape(2,3,4), W=np.arange(20,dtype=np.float64).reshape(4,5), bias=np.arange(5,dtype=np.float64), with no random draws. Predict Y=X@W+bias has shape (2,3,5), Y[1,2,4] approximately 1018, and Y.sum(axis=1,keepdims=True) shape (2,1,5). The standard product count is 2*3*5*(2*4-1)=210 scalar operations, plus 30 bias additions, with 25 shared parameters. Then use A=diag(1,1e-6), b=(1,0), db=(0,1e-8), all float64; x=solve(A,b), xp=solve(A,b+db). Measure \(\kappa_2\) approximately 1e6, relative solution change norm(xp-x,2)/norm(x,2) approximately 0.01, and original-system relative residual norm(b-A@xp,2)/norm(b,2) approximately 1e-8. These attain the bound by design; use absolute tolerance 1e-12 for the two ratios and entry, and relative tolerance 1e-12 for \(\kappa\). This changes the right-hand side, not the solver.Traps
- A successful broadcast can be wrong: subtracting a
(3,)target from a(3,1)prediction produces(3,3), not three residuals. Assert shapes before arithmetic. - Array
.Treverses all axes. For batched matrix transpose useswapaxes(-1,-2); transpose only the two matrix axes. - A small residual is backward-error evidence, not a guarantee of accurate \(x\). Conditioning belongs to the mathematical problem; stability belongs to the algorithm.
Practice
VMLS Exercise 10.11: trace, shapes and computation cost
Chapter 10, Exercise 10.11(a) to (d), printed p. 192
Derive the direct elementwise trace computation and its cost instead of forming the whole Gram product. The book's matrix norm here is Frobenius, not the induced 2-norm. The section mini-project should add batched shape assertions and a table of residual versus solution error as conditioning varies.
FNC Exercise 2.8.7: condition number of a product
Python edition, section 2.8.3, Exercise 2.8.7(a) and (b)
Use submultiplicativity and \((AB)^{-1}=B^{-1}A^{-1}\), then construct a strict-inequality example. State the induced norm and nonsingularity assumptions.
Self-test
-
What shapes result from
(2,1,3,4)@(1,5,4,6), and from swapping only the last two axes of the result?The matrix product is(2,5,3,6): 4 contracts and batch shapes(2,1)and(1,5)broadcast to(2,5). Swapping the matrix axes gives(2,5,6,3). Reversing every axis with.Twould instead give(6,3,5,2). -
Why does forming \(A^\top A\) square the 2-norm condition number, and which rank assumption is needed?
The SVD gives \(A^\top A=V\operatorname{diag}(\sigma_i^2)V^\top\). With full column rank all \(n\) eigenvalues are positive, so their extreme ratio is \((\sigma_{\max}/\sigma_{\min})^2\). A wide full-row-rank matrix still has a singular Gram matrix, so the finite full-column-rank formula cannot be used there. -
Where is invertibility used in the perturbation proof, and why does a residual give the same bound?
Invertibility lets \(A\Delta x=\Delta b\) become \(\Delta x=A^{-1}\Delta b\). Induced-norm inequalities bound the two products, and \(b\) nonzero permits relative normalization. If \(r=b-A\hat{x}\), then \(\hat{x}\) solves the fixed-\(A\) system with right-hand side \(b-r\); substitute \(\Delta b=-r\). No claim is made about the algorithm producing \(\hat{x}\). -
Why is reusing factors cheaper for several right-hand sides?
For dense \(n\) by \(n\) \(A\), one LU costs about \((2/3)n^3\). Each subsequent right-hand side requires two triangular solves, about \(2n^2\) operations. Recomputing the factorization for every vector repeats cubic work. Pivot comparisons and data movement are not captured by these arithmetic counts.
Completion criterion
Stop when
You predict the batched output and reduction shapes, derive the operation counts and both conditioning bounds, and explain how the worked case has a tiny residual but a one-percent solution change.
Mark passed, then open 2.1 Only mark this when you can do it without notes.Leads to
- 1.4 Branches, broadcasting and saved values in backward Apply the forward shape rules when deriving broadcast gradients.
- 2.1 Least squares: residuals, normal equations and stable solvers Apply the SVD proof of condition-number squaring to solver comparisons.
- 8.4 Precision formats, autocast and gradient scaling Use conditioning and arithmetic precision to reason about mixed-precision risk.
Comes after
- 1.7 SVD, pseudoinverse and low-rank approximationRead amplification and condition numbers from the singular spectrum.
Resources for this unit
- MIT 18.065 - Applied matrix bridge Use lectures 3-10 after 18.06 for SVD, PCA, low-rank approximation, matrix norms and solver difficulty.
- Introduction to Applied Linear Algebra Sample Chapters 1-13 and the extra exercises for computational practice in fitting, least squares and pseudoinverses.
- NumPy broadcasting and matmul shape rules Read General broadcasting rules and Broadcastable arrays, then
matmulNotes for stacked matrices and one-dimensional promotion. These close the pool's array-semantics gap. - Fundamentals of Numerical Computation: conditioning and operation counts Python edition, section 2.8 through Residual and backward error for the fixed-matrix perturbation bound; section 2.5.3 for triangular-solve and LU flop counts. Use Exercise 2.8.7 for a condition-number proof, not the unseeded Hilbert experiment.
Also in this section
- MIT 18.06 - Linear Algebra Follow lectures 1-10, 14-17, 21-22, 25, 27, 29-31 and 33 with matching problems. This is the ordered first pass.
- MIT 6.042J - Proof strand Study textbook Chapters 1 and 3, §§4.1-4.4 and §§5.1-5.4 with matching Unit 1 problems. Leave the rest of the course.
- Mathematics for Machine Learning Use Chapters 2-4 and §10.2 to check notation or re-derive a result. Do not treat it as a second course.
- Essence of Linear Algebra Open only when a geometric picture is missing; it does not replace derivation or problem solving.