Update for QR factorization least squares

I found after some research that the most numerically stable way to solve the least squares problem is through QR factorization. For $$n$$ number of observations and $$p$$ number of parameters it takes the following form:

$$Q_{n,p}R_{p,p} = X_{n, p}$$

Where $$R$$ is a square upper-triangular and $$Q$$ is orthogonal. And the solution is:

$$\mathbf{\hat{\beta}} = \left(\mathbf{X}^T \mathbf{X}\right)^{-1} \mathbf{X}^T z$$ $$\mathbf{\hat{\beta}} = R^{-1} Q^T z$$

I am interested on continuouly update my least squares model $$\hat{\beta}$$ as new data arrives $$x_{n+1},x_{n+2},x_{n+3}... \in R^{p}$$, we can define $$X^{*} = \begin{bmatrix} X \\ x_{n+1} \\ x_{n+2} \\ ... \end{bmatrix}$$ with added rows and items on $$y$$.

I found that is somehow called rank k matrix update or related to symmetric rank k update (SYRK) level 3 BLAS. I could not find more information on how formulate this on my LS QR factorization.

Any ideas or hints are really welcome.

SYRK is not really relevant here I think; it is just something else that happens to have the same name "rank k update".

In your case, you need to know how to update a QR factorization by inserting rows; a good reference is Golub, Van Loan, section 6.5.3: Appending or Deleting a Row.

Many computational environments have it already implemented for you, see e.g., Matlab's qrinsert, Python's scipy.linalg.qr_insert, Julia's QRupdate.jl, or this Fortran package.

• @thank you very much Federico – eusoubrasileiro Mar 2 '20 at 18:48
• unfortunately I guess this is not usually GPU optimized – eusoubrasileiro Mar 3 '20 at 11:11
• No, not that I know. I don't even know if there is anything optimized for adding multiple rows at the same time. – Federico Poloni Mar 3 '20 at 11:16

If you will be adding lots of rows, then you will want to use the tall skinny QR algorithm (TSQR) of Demmel et al, 2008,

https://arxiv.org/abs/0806.2159

This algorithm can be combined with the level 3 BLAS QR algorithm of Elmroth and Gustavson, 2000 to efficiently update the factorization with each new block of rows. It is best to save the new rows into a large block and do one update, rather than update with one row at a time.

A C implementation can be found in the GSL library,

https://www.gnu.org/software/gsl/doc/html/lls.html#tall-skinny-qr-tsqr-approach