qtyr#

Purpose#

Computes the orthogonal-triangular (QR) decomposition of a matrix \(X\) and returns \(Q'Y\) and \(R\).

Format#

{ qty, r } = qtyr(y, X)#

Parameters:

y (NxL matrix) – data
X (NxP matrix) – data

Returns:

qty (NxL matrix) – unitary matrix
r (KxP matrix) – upper triangular matrix. \(K = min(N,P)\).

Examples#

The QR algorithm is the numerically superior method for the solution of least squares problems:

loadm x, y;
{ qty, r } = qtyr(y, x);
q1ty = qty[1:rows(r), .];
q2ty = qty[rows(r)+1:rows(qty), .];

// LS coefficients
b = qrsol(q1ty, r);

// Residual sums of squares
s2 = sumc(q2ty^2);

Remarks#

Given \(X\), there is an orthogonal matrix \(Q\) such that \(Q'X\) is zero below its diagonal, i.e.,

\[\begin{split}Q′X = \begin{bmatrix} R \\ 0 \end{bmatrix}\end{split}\]

where \(R\) is upper triangular. If we partition

\[Q⁢ = \begin{bmatrix} Q_1 & Q_2 \end{bmatrix}\]

where \(Q_1\) has \(P\) columns, then

\[X = Q_1R\]

is the QR decomposition of \(X\). If \(X\) has linearly independent columns, \(R\) is also the Cholesky factorization of the moment matrix of \(X\), i.e., of \(X'X\). For most problems \(Q\) or \(Q_1\) is not what is required. Rather, we require \(Q'Y\) or \(Q_1'Y\) where \(Y\) is an NxL matrix (if either \(QY\) or \(Q_1Y\) are required, see qyr()). Since \(Q\) can be a very large matrix, qtyr() has been provided for the calculation of \(Q'Y\) which will be a much smaller matrix. \(Q_1'Y\) will be a submatrix of \(Q'Y\). In particular,

\[G = Q_1'Y = \text{qty}[1:P, .]\]

and \(Q_2'Y\) is the remaining submatrix:

\[H⁢ = Q_2'Y = \text{qty}[P+1:N, .]\]

Suppose that \(X\) is an NxK dataset of independent variables, and \(Y\) is an Nx1 vector of dependent variables. Then it can be shown that

\[b = R^{-1}G\]

and

\[\begin{split}s_j= \sum_{i=1}^{N−P}⁢H_{i,j}\\ ⁢j = 1,2,...L\end{split}\]

where b is a PxL matrix of least squares coefficients and s is a 1xL vector of residual sums of squares. Rather than invert \(R\) directly, however, it is better to apply qrsol() to

\[Rb⁢= Q_1′Y\]

For rank deficient least squares problems, see qtyre() and qtyrep().

Source#

qtyr.src