As a segue to my [prior post on this topic](http://stats.stackexchange.com/a/206345/67822) I want to share some tentative (albeit incomplete and flawed) exploration of the functions behind the linear algebra and related R functions. This is supposed to be a work in progress.

Part of the opaqueness of the functions has to do with the "compact" form of the Householder $\mathbf {QR}$ decomposition. The idea behind the Householder decomposition is to reflect vectors across a hyperplane determined by a unit-vector $\mathbf u$ as in the diagram below, but picking this plane in a purposeful way so as to project every column vector of the original matrix $\bf A$ onto the $\bf e_1$ standard unit vector.  The normalized norm-2 $1$ vector $\bf u$ can be used to compute the different Householder transformations $\mathbf{ I - 2\,  uu^T x}$.

[![enter image description here][1]][1]


The resultant projection can be expressed as

$\text{sign}(x_i=x_1)\times \lVert x \rVert \begin{bmatrix}1\\0\\0\\\vdots\\0\end{bmatrix}+\begin{bmatrix}x_1\\x_2\\x_3\\\vdots\\x_m\end{bmatrix}$

The vector $\bf v$ represents the difference between the column vectors $\bf x$ in the matrix $\bf A$ that we want to decompose and the vectors $\bf y$ corresponding to the reflection across the subspace or "mirror" determined by $\bf u$.

The method used by LAPACK liberates the need for storage of the first entry in the Householder reflectors by turning them into $1$'s. Instead of normalizing the vector $\bf v$ to $\bf u$ with $\lVert u\rVert= 1$, it is just the fist entry that is converted to a $1$; yet, these new vectors - call them $\bf w$ can still be used as a directional vectors. 

The beauty of the method is that given that $\bf R$  in a $\bf QR$ decomposition is upper triangular, we can actually take advantage of the $0$ elements in $\bf R$ below the diagonal to fill them in with these $\bf w$ reflectors. Thankfully, the leading entries in these vectors all equal $1$, preventing a problem in the "disputed" diagonal of the matrix: knowing that they are all $1$ they don't need to be included, and can yield the diagonal to the entries of $\bf R$.

The "compact QR" matrix in the function `qr()$qr` can be understood as roughly the addition of the $\mathbf R$ matrix and the lower triangular "storage" matrix for the "modified" reflectors. 

The Householder projection will still have the form $\mathbf{ I - 2\,  uu^T x}$, but we won't be working with $\bf u$ ($\lVert \bf x \rVert=1$), but rather with a vector $\bf w$, of which only the first entry is guanteed to be $1$, and $\mathbf{ I - 2\,  \frac{w}{\lVert w \rVert} \frac{w^T}{\lVert w \rVert} x}$

$\large u[2:\text{length(u)}]/\tau$ 

with 

$\large u = \frac{\text{sign}(x_i=x_1)\times \lVert x \rVert \begin{bmatrix}1\\0\\0\\\vdots\\0\end{bmatrix}+\begin{bmatrix}x_1\\x_2\\x_3\\\vdots\\x_m\end{bmatrix}}{\text{sign}(x_i=x_1)\times \lVert x \rVert + x[1]}$ and $\large \tau = \frac{1 + u[2:\text{length(u)}]^T\,u[2:\text{length(u)}]}{2}$.

    $PseudoQ
           [,1]   [,2]  [,3] [,4]
    [1,]  1.000  0.000 0.000    0
    [2,] -0.148  1.000 0.000    0
    [3,]  0.939 -0.675 1.000    0
    [4,]  0.099  0.042 0.623    1

In this way $tau$ or other derived scalars, such as $\rho$ (in linked code page) can be easily stored to reproduce $Q$, while in the lower $0$ part of $R$ these "decapitated" $u$ vectors are neatly stored (in red below):


[![enter image description here][2]][2]

The equivalence between the manual and R functions are self-explanatory in the link to the code attached in Github.

Although these computations are meant to facilitate working with block matrices, the intuition can be gained with a simple random matrix. All the code that follows is [here](https://github.com/RInterested/SIMULATIONS_and_PROOFS/blob/master/HOUSEHOLDER%20QR).

  [1]: https://i.sstatic.net/eOHkX.png
  [2]: https://i.sstatic.net/LmydM.png