[SciPy-dev] Why does orth use svd instead of QR ?
Fri Feb 5 02:29:51 CST 2010
Charles R Harris wrote:
> On Fri, Feb 5, 2010 at 12:47 AM, Charles R Harris
> <email@example.com <mailto:firstname.lastname@example.org>> wrote:
> On Fri, Feb 5, 2010 at 12:18 AM, David Cournapeau
> <email@example.com <mailto:firstname.lastname@example.org>> wrote:
> Charles R Harris wrote:
> > On Thu, Feb 4, 2010 at 8:45 PM, David Cournapeau
> <email@example.com <mailto:firstname.lastname@example.org>
> > <mailto:email@example.com
> <mailto:firstname.lastname@example.org>>> wrote:
> > Hi,
> > I wanted to know if there was a rationale for using svd to
> > orthonormalize the columns of a matrix (in scipy.linalg).
> > methods are likely to be much faster, and I thought this
> was the
> > standard, numerically-stable method to orthonormalize a
> basis ? If the
> > reason is to deal with rank-deficient matrices, maybe we
> could add an
> > option to choose between them ?
> > QR with column rotation would deal with rank-deficient
> matrices and
> > routines for that are available in LAPACK
> > <http://netlib.org/lapack/lug/node42.html>. The SVD was
> probably used
> > because it was available. The diagonal elements of the R
> matrix can
> > somewhat take the place of the singular values when column
> rotation is used.
> So would be it ok to use this column-rotated QR in place of svd for
> every case in orth ? I would have to check that QR with column
> is still significantly faster than svd, but I would surprised if
> if were
> not the case. QR has also the advantage of being implemented in
> already contrary to eigen/svd solvers,
> I don't know how the two methods compare in practice. SVD algorithms
> generally use iterated QR reductions in their implementation, so QR
> reductions can't be worse numerically. But the SVD probably provides
> a better metric for rank determination. A google search turns up
> some literature on the subject that I can't access from home.
> OK, here's a good reference
> <http://www.math.sjsu.edu/%7Efoster/rank/rank_revealing_s.pdf>. A quick
> look seems to indicate that the SVD is the way to go.
AFAIK, SVD is indeed the way to go, but do we really need this for the
orth function ? I am wrong to think that orthonormalizing a matrix of
linearly independent vectors is the most common usage for orth ? The
difference in terms of speed is really significant (for example, svd of
a 2000x100 matrix takes ~1.9 second vs 0.1 s for QR).
More information about the SciPy-Dev