3.4 Practice Problems
As in 2.7, predict before you run. Several of these ask you to verify a property rather than compute a number—that habit is the one worth building, because it is how you catch a wrong answer that looks right.
Problems 1–9 are open-ended: write code in the cell, and judge your own answer against the property you were asked to check. Problems 10–12 are auto-graded and look ahead to specific later lectures.
Solving systems
Solve A\mathbf{x} = \mathbf{b} two ways, with solve and with inv. Confirm they agree, then time both on a 500\times 500 system and say which you would use in practice.
This system has no solution. Predict what solve does, then what lstsq does, and explain the difference in terms of the column space.
Projection
Build the projection matrix P onto the column space of A. Verify all three properties without looking them up: P is symmetric, P^2 = P, and P\mathbf{v} = \mathbf{v} for any \mathbf{v} already in the column space.
Fit a straight line to the points below by least squares. Then show that the residual vector is orthogonal to both columns of the design matrix, and explain why that was guaranteed before you computed it.
The trace of a projection matrix equals the dimension of the space it projects onto. Check this for the P from Problem 3, then predict the trace for a design matrix with 5 linearly independent columns.
Decompositions
Compute A = QR. Confirm Q^\top Q = I and QR = A, then compute the least-squares coefficients from Q and R alone, without forming A^\top A.
Take a symmetric matrix and compute its eigendecomposition with both eig and eigh. Compare the outputs and say why eigh is the right choice here.
Truncate the SVD of X at each rank r = 1,\dots,6 and record \|X - \tilde{X}\|_F. Plot the six errors against r, and on the same axes plot \sqrt{\sum_{k>r}\sigma_k^2}. Say what you expect the plot to look like before you draw it.
Where this is going
Centre X by column, form the covariance matrix, and compute its eigenvalues with eigh. Separately, compute the singular values of the centred X with svd. The two are related by a simple formula—find it empirically, then say what it means for PCA.
Two predictors are nearly identical. Compute the condition number of A^\top A, then fit least squares and look at the coefficients. Perturb b by 0.001 in one entry and refit. Report how much the coefficients moved, and how much the fitted values moved.
This is multicollinearity, and it is the problem ridge regression exists to solve on Nov 3.
Auto-graded
These three are graded, and each one is a preview of a specific later lecture. Predict the answer, then check it.
pair_dists(Y) returns all pairwise Euclidean distances among the rows of Y. For a cloud of 200 uniform points, compute the coefficient of variation of those distances—the standard deviation divided by the mean—in 2 dimensions and in 100 dimensions, and report the ratio \mathrm{cv}_{100}/\mathrm{cv}_{2}.
The coefficient of variation is the right summary here because it is scale-free: distances get bigger with dimension no matter what, and the question is whether they get relatively more alike.
answer = (d100.std() / d100.mean()) / (d2.std() / d2.mean())
The mean distance grows like \sqrt{d}, but the spread grows much more slowly, so the cloud looks more and more like a thin shell at a fixed radius. The remedy on Nov 17 is not a better distance function; it is fewer dimensions—which is what PCA provides.
A is the near-collinear design matrix from Problem 9 and S holds its singular values, which are (7.96,\ 0.795,\ 1.45\times10^{-4}). Ordinary least squares weights the k^{\text{th}} singular direction by 1/\sigma_k, i.e. by (0.126,\ 1.259,\ 6901)—look at that third number. Ridge with penalty \lambda replaces the weight by the filter factor \sigma_k/(\sigma_k^2+\lambda). Compute all three filter factors for \lambda = 0.01.
answer = S / (S**2 + lam)
When \sigma_k^2 \gg \lambda the factor is \approx 1/\sigma_k, unchanged from OLS. When \sigma_k^2 \ll \lambda it is \approx \sigma_k/\lambda \approx 0. The penalty \lambda is therefore a threshold on \sigma^2, separating directions the data determined from directions it did not. This is the sentence you will want on Nov 3, and it is only sayable because of the SVD.
X holds 300 points in \mathbb{R}^{50}: three genuine clusters living in the first two coordinates, plus 48 coordinates of pure noise. sep(Y) measures how separated the clusters are—mean between-cluster distance divided by mean within-cluster distance—so bigger is better and 1 means the clusters are invisible.
On the raw 50-dimensional data, sep(X) is only 1.17: the clusters are all but undetectable. Centre X, project it onto its top two principal components, and compute sep of the result.
Z has shape (300, 2): the coordinates of each observation in the principal-component basis. Note that the top two components carry only about 21\% of the total variance, so an “explained variance” rule of thumb alone would not have told you to stop at two—the structure is in those components even though most of the variance, spread thinly across 48 noise directions, is not.