The regression model: a target plus noise

Until now we only had a data table and looked for structure in it. Now each row also has an answer, and we want a rule that predicts the answer for new rows.

\(X_1,\dots,X_d\)

predictors (inputs), \(d\) of them

\(Y\)

response (the thing to predict)

\(f(\mathbf X)\)

the part of \(Y\) explained by the predictors

\(\varepsilon\)

random error, independent of \(\mathbf X\): measurement noise plus everything not measured (latent variables)

\(n\)

number of training rows

Delivery time depends on distance, but also on traffic, weather and the driver. What we do not measure ends up in \(\varepsilon\). The slides' examples: food delivery (distance, stops → time), power grid (temperature, hour → load), real estate (area, rooms → price), language models (compute → test loss).

The linear model: bias plus weights

We guess that the answer is a weighted sum of the inputs plus a constant. With one input it is a straight line, with two it is a flat plane.

\(b\)

bias or intercept: the prediction when every input is 0

\(w_i\)

regression coefficient (weight) of \(X_i\): how much \(\hat y\) changes per unit of \(X_i\)

\(\hat y_i\)

predicted value for row \(i\) (hat = fitted, not observed)

\(\epsilon_i\)

residual: observed minus predicted, for row \(i\)

Least squares: pick the line with the smallest squared error

There are many possible lines. We choose the one whose total squared miss, added over all points, is smallest.

Because SSE is a smooth bowl in \(b\) and \(\mathbf w\), setting its derivatives to zero gives the minimum directly. Every formula on this page comes from doing exactly that.

One predictor: intercept and slope

For one input, the best line always goes through the point of averages, and its slope is the covariance of the two variables divided by the variance of the input.

\(\mu_X, \mu_Y\)

sample means of the input and the response

\(\sigma_{XY}\)

sample covariance of \(X\) and \(Y\)

\(\sigma_X^2\)

sample variance of \(X\)

Recipe: fit a line by hand

1. Compute \(\mu_X\) and \(\mu_Y\).
2. Make a table with columns \(x_i - \mu_X\), \(y_i - \mu_Y\), \((x_i - \mu_X)^2\), \((x_i-\mu_X)(y_i-\mu_Y)\).
3. \[ w = \frac{\sum (x_i - \mu_X)(y_i - \mu_Y)}{\sum (x_i - \mu_X)^2} \]
4. \[ b = \mu_Y - w\,\mu_X \]
5. Write the line \(\hat y = b + w x\) and use it to predict.

Delivery times: fit the line by hand and predict 7 km

A delivery company records five orders in Doha. How long should a 7 km order take?

Order	\(x\) (km)	\(y\) (min)
1	2	15
2	4	16
3	6	24
4	8	30
5	10	35

1. \[ \mu_X = \frac{30}{5} = 6 \]
2. \[ \mu_Y = \frac{120}{5} = 24 \]
3. Deviation table:

   \(x_i\)	\(y_i\)	\(x_i - \mu_X\)	\(y_i - \mu_Y\)	\((x_i-\mu_X)^2\)	\((x_i-\mu_X)(y_i-\mu_Y)\)
   2	15	−4	−9	16	36
   4	16	−2	−8	4	16
   6	24	0	0	0	0
   8	30	2	6	4	12
   10	35	4	11	16	44
   			sum	40	108

4. \[ w = \frac{108}{40} = 2.7 \]
5. \[ b = 24 - 2.7 \times 6 = 7.8 \]
6. \[ \hat y = 7.8 + 2.7\,x \]
7. \[ \hat y(7) = 7.8 + 2.7 \times 7 = 26.7 \text{ min} \]

Interpretation (s11): each extra km adds about 2.7 minutes; 7.8 min is the fixed time per order (pickup, parking). Predicting at 30 km is extrapolation: there is no data there.

Residuals and their two zero-sum identities

A residual is how far each point sits above or below the line. For a least-squares line, the residuals always add up to zero and are not correlated with the input.

Setting \(\partial\,\mathrm{SSE}/\partial b = 0\) is the same as \(\sum \epsilon_i = 0\). Setting \(\partial\,\mathrm{SSE}/\partial w = 0\) is the same as \(\sum x_i \epsilon_i = 0\). The geometry card residuals are orthogonal says the same thing with vectors.

Root mean squared error: a typical miss in the units of \(Y\)

SSE grows with the number of points and is in squared units. RMSE fixes both, so it tells you how far off a typical prediction is.

Residual plots: what a good fit looks like

Plot each residual against its prediction. If you can see a shape, the straight line is missing something, no matter how small the error number is.

Iris: petal width from petal length

The same recipe on 150 flowers. The numbers are bigger, but the only inputs are the two means, one variance and one covariance.

Two views of the data: rows as points, columns as vectors

So far each order was a dot on a 2-D plot. Now flip it: each variable becomes one long vector with one entry per order, and fitting a line means choosing one vector.

\(\mathbf 1\)

the vector of \(n\) ones (it carries the intercept)

\(X, Y\)

capital letters now mean the whole column of \(n\) values

\(\operatorname{span}\{\mathbf 1, X\}\)

every vector of the form \(b\mathbf 1 + wX\): a 2-D plane inside \(\mathbb R^n\)

Least squares = orthogonal projection of \(Y\) onto the plane

The data vector is not in the plane because the data are noisy. The closest point in the plane is found by dropping a perpendicular, and that closest point is the fitted vector.

Projecting onto a unit vector was card “Projection” on the Lectures 3–4 page. The only new thing is that the target is a plane with two directions, and they are not at right angles, so we first split them.

Recipe: the “geometric approach” a question may ask for

1. Write \(\mathbf 1\), \(X\) and \(Y\) as column vectors of length \(n\).
2. Form \(A = [\mathbf 1 \; X]\) and compute \(A^T A\) and \(A^T Y\).
3. Solve \((A^T A)\,(b, w)^T = A^T Y\) with the 2×2 inverse (next part has the formula).
4. \[ \hat Y = A\,(b, w)^T = b\,\mathbf 1 + w\,X \]
5. Check: \(\boldsymbol\epsilon = Y - \hat Y\) sums to zero and \(X^T \boldsymbol\epsilon = 0\).

Absorb the bias: one weight vector for all \(d + 1\) parameters

Add a column of ones to the data, and the intercept becomes just another weight. Then one matrix formula covers one predictor, two, or a hundred.

\(\tilde D\)

the data matrix with a column of ones in front (tilde = augmented)

\(\tilde{\mathbf w}\)

bias first, then the \(d\) weights

\(\tilde{\mathbf x}_i\)

row \(i\) of \(\tilde D\), i.e. \((1, x_{i1}, \dots, x_{id})^T\)

With \(d = 1\), \(\tilde D\) is exactly the matrix \(A = [\mathbf 1 \; X]\) from the projection card. HW2 Q6 calls it \(A\).

Normal equations: the closed-form solution

Write SSE with matrices, set its gradient to zero, and the best weights come out of one matrix inverse. This is the formula the exam will give you.

\(\tilde D^T \tilde D\)

\((d+1) \times (d+1)\) matrix of sums of products of the columns

\(\tilde D^T Y\)

vector of sums of products of each column with \(Y\)

Recipe: normal equations by hand

1. Build \(\tilde D\) by putting a column of ones before the predictors.
2. Compute \(\tilde D^T \tilde D\): entry \((j,k)\) is the dot product of columns \(j\) and \(k\). The top-left entry is \(n\).
3. Compute \(\tilde D^T Y\): entry \(j\) is the dot product of column \(j\) with \(Y\).
4. Invert \(\tilde D^T \tilde D\) (2×2 by formula; 3×3 only if the question makes it easy) and multiply.
5. Read off \(b\) first, then \(w_1, \dots, w_d\).

Hand calculation with centred data

Subtract the means first and the intercept drops out of the matrix. With two predictors you only need to invert a 2×2 matrix, then get the bias from the means.

\(\bar D\)

\(n \times d\) matrix of centred predictors (no column of ones)

\(\boldsymbol\mu_X\)

vector of the \(d\) predictor means

\(\Sigma\)

covariance matrix of the predictors, with \(1/n\)

Recipe: two predictors by hand

1. Compute the means \(\mu_1, \mu_2, \mu_Y\) and centre all three columns.
2. \[ \bar D^T \bar D = \begin{bmatrix} \sum \bar x_{i1}^2 & \sum \bar x_{i1}\bar x_{i2} \\ \sum \bar x_{i1}\bar x_{i2} & \sum \bar x_{i2}^2 \end{bmatrix} \]
3. \[ \bar D^T \bar Y = \begin{bmatrix} \sum \bar x_{i1} \bar y_i \\ \sum \bar x_{i2} \bar y_i \end{bmatrix} \]
4. Invert the 2×2 with the formula and multiply to get \(\mathbf w = (w_1, w_2)^T\).
5. \[ b = \mu_Y - w_1 \mu_1 - w_2 \mu_2 \]

Iris with two predictors: a fitted plane

With sepal length added as a second input, the model is a plane and the error barely improves. The numbers show what the normal equations look like at full size.

Top-left is \(n = 150\). The first row holds the column sums \(\sum x_{i1} = 876.5\) and \(\sum x_{i2} = 563.8\). The rest are sums of squares and cross products. You will not be asked to invert a 3×3 like this by hand; the centred route reduces it to 2×2.

Correlated predictors: why a weight can flip sign

When two inputs carry almost the same information, the model cannot tell which one deserves the credit, so individual weights become unstable and hard to interpret.

High correlation makes \(\bar D^T \bar D\) nearly singular: its determinant \(ae - c^2\) is small, so the inverse has large entries and small changes in the data swing the weights a lot.

In software, and a warning about training error

Software does the whole fit in three lines. But a small error on the training data does not prove the model will predict well on new data.

from sklearn.datasets import load_iris
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_squared_error

X = load_iris().data
D, y = X[:, [0, 2]], X[:, 3]        # sepal & petal length -> petal width
model = LinearRegression().fit(D, y)
print(model.intercept_, model.coef_)   # b, (w1, w2)
y_hat = model.predict(D)
print(mean_squared_error(y, y_hat) ** 0.5)     # RMSE

The old model is a special case of the new one (set the new weight to 0). Least squares picks the best of a bigger set, so it can only do at least as well. That is also why training error is a poor judge.

Which formula, when