[R] Is there a bug in qr(..,LAPACK=T)

[Ravi Varadhan]

As the help page for "qr" says, LAPACK does not attempt to detect 
linear dependencies or rank deficiencies, so you should not use the 
value of "rank" obtained with argument, LAPACK = TRUE. Computing the 
rank of a matrix using finite precision is difficult, as the example on 
the help page for "qr" shows using Hilbert matrix order 9.  The rank 
can change depending on "tolerance" option, which is actually not used 
if LAPACK = TRUE.

Ravi.

----- Original Message -----
From: Mike Meyer <mikem at salter-point.com>
Date: Wednesday, July 16, 2003 12:54 pm
Subject: [R] Is there a bug in qr(..,LAPACK=T)

> The following snippet suggests that there is either a bug in 
> qr(,LAPACK=T), or some bug in my understanding.   Note that the 
> detected rank is correct (= 2) using the default LINPACK qr, but 
> incorrect (=3) using LAPACK.   This is running on Linux Redhat 
> 9.0, using the lapack library that comes with the Redhat 
> distribution.   I'm running R 1.7.1 compiled from the source.  If 
> the bug is in my understanding (or in the Redhat 9.0 libraries or 
> compiler) I would much appreciate some enlightenment.
> Thanks, --Mike
> 
> > X
>     [,1] [,2] [,3]
> [1,]    1    1    1
> [2,]    1    2    1
> [3,]    1    3    1
> [4,]    1    4    1
> > qr(X)
> $qr
>     [,1]       [,2] [,3]
> [1,] -2.0 -5.0000000   -2
> [2,]  0.5 -2.2360680    0
> [3,]  0.5  0.4472136    0
> [4,]  0.5  0.8944272    0
> 
> $rank
> [1] 2
> 
> $qraux
> [1] 1.5 1.0 0.0
> 
> $pivot
> [1] 1 2 3
> 
> attr(,"class")
> [1] "qr"
> 
> > qr(X,LAPACK=T)
> $qr
>           [,1]       [,2]          [,3]
> [1,] -5.4772256 -1.8257419 -1.825742e+00
> [2,]  0.3087742 -0.8164966 -8.164966e-01
> [3,]  0.4631613 -0.3270981 -1.378276e-16
> [4,]  0.6175484 -0.7892454  9.055216e-01
> 
> $rank
> [1] 3
> 
> $qraux
> [1] 1.182574 1.156135 1.098920
> 
> $pivot
> [1] 2 1 3
> 
> attr(,"useLAPACK")
> [1] TRUE
> attr(,"class")
> [1] "qr"
> >
> 
> -- 
> 
> Mike Meyer,  Seattle WA
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://www.stat.math.ethz.ch/mailman/listinfo/r-help
>