Algorithm for minimization of sum of squares in regression packages

I'd like to know what technique or techniques are used by regression packages, (in particular the lm function of R) to minimize the sum of squares.

Thanking in anticipation.

-

No, it doesn't use gradient descent to fit linear models.

Linear least squares has an explicit solution.

If we ignore weights, and the possibility of multiple $y$'s, and just deal with "plain" multiple regression:

$E(y) = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \ldots + \beta_p x_p + \epsilon$

$\quad\quad\,\,=X\beta+\epsilon$

Then attempting to find the argmin of the sum of squares of errors lead to the least squares normal equations, $(X^TX)\hat\beta=X^Ty$ which have the algebraic solution $\hat\beta=(X^TX)^{-1}X^Ty$ -- but lm doesn't actually compute that.

What most regression programs do instead (lm included) is to compute the QR decomposition of $X$, and then the normal equations become:

$(R^TQ^TQR)\hat\beta=X^Ty$, but $Q^TQ=I$, so

$(R^TR)\hat\beta=(R^T Q^T)y$ which can be recast as

$R^T(R\hat\beta)=R^T (Q^Ty)$

And then (skimming over quite a few details*) the fact that $R$ is upper triangular is exploited to solve that system efficiently.

* (including the use of pivots/permutation matrices, the big-R/little-R dichotomy, simplifying the above further before solving, and a bunch of other issues)

If you search on QR decomposition least squares you should find sets of notes that lay out the full details (but you'll likely have to learn a number of things before it's all clear).

A classic reference is Golub and Van Loan's Matrix Computations.

While this is - more or less - the way most least squares regression code works these days, you may find some that use either Choleski decomposition or singular value decomposition of the $X^TX$ matrix, or in a few cases some other algorithm.

-