Linear and logistic regression from scratch

One model, three ways to fit it

A linear model predicts ŷ = Xw + b. Fitting it means choosing w and b to minimise a loss, and there are three standard routes. Interviews use all three, often in the same hour, because each exposes a different kind of misunderstanding.

• Closed form. For squared error the minimum solves a linear system. Exact, and the right answer when there are a few thousand features or fewer.
• Gradient descent. Works for any differentiable loss, scales to data that does not fit in memory, and is the only route open to logistic regression.
• A library call. np.linalg.lstsq or scikit-learn. Fine in practice, and not what the interviewer is asking for.

The normal equations

For the loss ‖Xw − y‖², the gradient is 2Xᵀ(Xw − y). Setting it to zero gives XᵀX w = Xᵀy. Adding a ridge penalty λ‖w‖² adds 2λw to the gradient, so

(XᵀX + λI) w = Xᵀy.

For λ > 0 the matrix on the left is symmetric positive definite, so the system always has exactly one solution, even when there are more features than observations. That is the practical reason ridge exists, before any argument about overfitting.

Solve the system; do not invert the matrix. np.linalg.solve factorises once and is both faster and more accurate than np.linalg.inv(A) @ b, and an interviewer who sees inv will ask why.

Gradient descent, vectorised

With the mean squared error L = (1/n)‖Xw + b − y‖² and residuals r = Xw + b − y,

∂L/∂w = (2/n) Xᵀr,   ∂L/∂b = (2/n) Σ r_i.

One step is three lines of NumPy, with no loop over examples. Shapes are the quickest check: X is (n, d), r is (n,), so Xᵀr is (d,), the same shape as w. If your gradient has the shape of the data rather than the parameters, it is wrong.

r = X @ w + b - y
w = w - lr * (2 / n) * (X.T @ r)
b = b - lr * (2 / n) * r.sum()

Both gradients are computed from the same residuals, before either parameter moves. Updating w and then recomputing the residuals for b is a different algorithm, and it will not reproduce anyone else's numbers.

Logistic regression

Replace the identity with the sigmoid, p = σ(Xw + b), and squared error with binary cross-entropy. The gradient is the reason the pairing is used: the derivative of the loss with respect to the logit is simply p − y, so

∂L/∂w = (1/n) Xᵀ(p − y) + λw,   ∂L/∂b = mean(p − y)

for the mean loss with a penalty (λ/2)‖w‖². The penalty does not touch the bias; penalising it would pull every prediction towards one half for no reason. There is no closed form, which is why this is the standard first exercise in gradient descent.

What goes wrong in practice

• Learning rate. For least squares, gradient descent converges only if the learning rate is below 2 divided by the largest eigenvalue of (2/n)XᵀX. Above it the loss grows geometrically. Unscaled features make that eigenvalue enormous; standardise first.
• Separable data. If a line separates the classes perfectly, unpenalised logistic regression has no minimum: the weights grow for ever. A small λ fixes it.
• Shapes. A y of shape (n, 1) against predictions of shape (n,) broadcasts to an (n, n) matrix and gives a wrong answer silently. Flatten early.