Lecture 13: Randomized Matrix Multiplication

What this lecture covers

For very large matrices, computing the exact product A * B can be too expensive, so this lecture introduces a randomized shortcut: sample a few column-row pairs of A and B, weight and average them, and get an approximation that is correct on average. The lecture works through the statistics needed to make this rigorous, then derives the sampling probabilities that make the approximation as accurate as possible.

This is the first lecture in the course to lean on mean and variance, framing randomized linear algebra as a probability problem. By the end, you can explain why sampling with probability proportional to the squared norm of each column-row pair gives the correct expected product and the smallest possible variance, and you'll have seen Lagrange multipliers used to prove that claim.

Key ideas

• Randomized matrix multiplication: instead of forming the full product A * B, you sample individual rank-one pieces A_j * B_j^T (a column of A times the corresponding row of B) with chosen probabilities, then average the samples.
• Mean and variance: the mean is the probability-weighted average of the possible outcomes; the variance measures the average squared distance of outcomes from that mean, computed either as sum p_i * (output_i - mean)^2 or the equivalent sum p_i * output_i^2 - mean^2.
• Unbiased sampling: dividing each sample by s * p_j (where s is the number of samples) makes the expected value of the random approximation equal to the true product A * B, for any valid choice of probabilities.
• Norm-squared weighting: the probability of picking column-row pair j should be proportional to norm(A_j) * norm(B_j) (norm squared when B = A^T), so bigger columns are sampled more often.
• Variance as the quantity to minimize: getting the mean right is not enough; the goal is to choose probabilities that make the variance of the approximation as small as possible.
• Lagrange multipliers: since the probabilities must sum to one, the optimization is a constrained minimization, solved by adding a lambda term for the constraint and setting partial derivatives to zero.
• Result: the derivative computation confirms that the optimal probabilities are exactly the norm-squared weights introduced earlier, so that initial guess turns out to be provably best.

Walkthrough

Setting up the idea (0:01)

Strang introduces randomized linear algebra for very large matrices: instead of computing A * B directly, you sample columns of A and matching rows of B, weight them by chosen probabilities, and add up the samples. He previews the plan: decide probabilities, compute the mean of the process, show it equals the true product, then compute the variance and minimize it by choosing the best probabilities.

A small practice example (4:09)

To recall what mean and variance mean, Strang works a toy case: a 1-by-2 matrix, sampled twice with equal probability (1/2 and 1/2) for its two columns. He computes the mean of this process by hand and confirms it equals the correct product, then sets up the variance calculation for the same example, reviewing both standard formulas for variance along the way.

Computing the variance for the example (11:31)

Working through the arithmetic, Strang finds a nonzero variance for the toy example, since no single sample ever equals the exact answer. He uses this result to make a key observation: if one column is much larger than the other, sampling both equally wastes effort, and weighting the larger column more heavily would shrink the variance.

Norm-squared probabilities (21:58)

Strang generalizes: probabilities for choosing column-row pair j should be proportional to the norm of that column times the norm of the corresponding row (norm squared in the symmetric case). Since these raw values don't sum to one, he rescales by their total sum C. He also mentions an alternative strategy of mixing columns into similar sizes first and sampling uniformly, but chooses to develop the norm-weighted approach.

Mean and variance of the general sampling process (28:07)

With the sampling rule set up (divide each sample by s times its probability), Strang shows the mean of the randomized approximation comes out exactly to A * B, regardless of which probabilities are used. He then turns to the harder part: computing the variance of this general process, which involves substituting the chosen probabilities into the variance formula and simplifying.

Minimizing the variance with Lagrange multipliers (38:31)

Strang poses the final optimization: minimize the variance expression, a sum of (A_j B_j)^2 / p_j terms, subject to the probabilities summing to one. He applies Lagrange's method, adding a multiplier lambda for the constraint, then takes partial derivatives with respect to each p_j and with respect to lambda and sets them to zero.

Solving for the optimal probabilities (44:46)

Working through the derivative equations, Strang solves for p_j and shows it must be proportional to the norm-squared quantity A_j B_j, confirming that the weighting scheme introduced earlier is in fact the variance-minimizing choice. He closes by summarizing the overall procedure: compute column and row lengths first, use them to set sampling probabilities, and only then run the randomized sampling.

Before you watch

• Be comfortable with matrix multiplication as a sum of rank-one column-times-row pieces, covered in earlier lectures on the four ways to multiply matrices.
• Basic familiarity with mean and variance is helpful, though Strang reviews both from scratch in this lecture.
• Knowing the idea of Lagrange multipliers for constrained optimization will make the final derivation easier to follow, though it is also explained as it's used.

Check your understanding

1. Why does making the mean of the randomized approximation correct not guarantee a good approximation?
2. In the sampling scheme, why must each sampled term A_j B_j^T be divided by s * p_j rather than just p_j?
3. Why does weighting the sampling probability by the norm of the column and row reduce the variance compared to uniform sampling?
4. What role does the Lagrange multiplier lambda play in finding the optimal probabilities, and what equation do you get by differentiating with respect to it?
5. If one column of A is much larger than the others, how should its sampling probability compare to the rest, and why?