[Feature Request] QR Matrix Decomposition using Gram-Schmidt

[LibraChris]

Introduce an option to utilize the Gram-Schmidt process for QR matrix decomposition alongside the existing Householder transformation method. This would allow users to choose the decomposition method based on their specific requirements, providing flexibility and potentially improving the overall usability of the library.

The Gram-Schmidt process and the Householder transformation are both methods used for QR matrix decomposition, but they differ in their computational approach and resulting matrix dimensions:

1. Gram-Schmidt Process:

   • Input: An $( m \times n $) matrix, where $m$ is the number of rows and $n$ is the number of columns.
   • Output: Two matrices $Q$ and $R$, where $Q$ is an $m \times n$ orthogonal matrix (i.e., $Q^TQ = I$) and $R$ is an $n \times n$ upper triangular matrix.

2. Householder Transformation:

   • Input: An $m \times n$ matrix, where $m$ is the number of rows and $n$ is the number of columns.
   • Output: Two matrices $Q$ and $R$, where $Q$ is an $m \times m$ orthogonal matrix (i.e., $Q^TQ = I$) and $R$ is an $m \times n$ upper triangular matrix.

In summary, the main difference lies in the dimensions of the orthogonal matrix $Q$: Gram-Schmidt produces an $m \times n$ orthogonal matrix, while the Householder transformation yields an $m \times m$ orthogonal matrix.

[bvenn]

Summary

As you mentioned there are two result types of a QR decomposition: (i) full/complete form QR (Householder) and (ii) reduced form QR (Gram-Schmidt)^1,2,3. Gram-Schmidt seems to be efficient and easy to compute, but is numerically unstable^4,5. There seems to be a modified Gram-Schmidt method with increased accuracy⁶. This may be worth looking at in future. The reduced form seems to be unable to handle rank-deficient matrices.

Unfortunately you cannot convert the full QR into the reduced form (at least I was not able to find a reference that shows how to do it). If you want, you can open a PR adding the Gram-Schmidt QR decomposion. As it's less numerically stable than the existing version I would suggest to call it QR_GramSchmidt and add detailed XML-documentation to both versions explaining the differences :+1:

References:

1. https://de.wikipedia.org/wiki/QR-Zerlegung#Definition
2. https://qph.cf2.quoracdn.net/main-qimg-c19b21f42d093963b26bf5784c8aeafb
3. https://studyflix.de/mathematik/qr-zerlegung-1786 Section "Allgemeine QR Zerlegung"
4. https://en.wikipedia.org/wiki/QR_decomposition#Advantages_and_disadvantages
5. https://blogs.mathworks.com/cleve/2016/07/25/compare-gram-schmidt-and-householder-orthogonalization-algorithms
6. https://en.wikipedia.org/wiki/Gram%E2%80%93Schmidt_process#Numerical_stability

[kMutagene]

@LibraChris @bvenn

Depending on how hard and time consuming a managed implementation would be, I would suggest to first go the route of adding bindings to the respective LAPACK method and using that instead for starters. I went the same route with the thin SVD.

Having a managed version in F# has value, however there is a reason why these linear algebra routines are used under the hood by almost any high performance stats/ML library - they are highly optimized and tested. So depending on the use case (and how fast you need it) i would highly suggest going this route. If i can give more pointers/help on that feel free to ask. As a first pointer, it looks to me like this is the needed subroutine: https://www.netlib.org/lapack/explore-html/d0/da1/group__geqrf_gade26961283814bb4e62183d9133d8bf5.html#gade26961283814bb4e62183d9133d8bf5