RE: [R] rank of a matrix

From: Ravi Varadhan <rvaradha_at_jhsph.edu>
Date: Thu 05 May 2005 - 08:03:18 EST


For the example that Duncan just gave, using the "matrix.rank" function of Spencer Graves (which uses singular value decomposition) I obtained the following result:

> exponent <- -(7:16)
> eps <- 10^exponent
> sapply(eps,mat=h9times4,function(x,mat)matrix.rank(mat,x))
 [1] 6 7 7 8 8 9 9 9 9 9

This tells me that the correct rank should be 9, since the rank stabilizes for smaller tolerances. I realize that this may not work generally, and one could create counter-examples to defeat this strategy.

Best,
Ravi.



Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: rvaradhan@jhmi.edu

> -----Original Message-----
> From: r-help-bounces@stat.math.ethz.ch [mailto:r-help-
> bounces@stat.math.ethz.ch] On Behalf Of Duncan Murdoch
> Sent: Wednesday, May 04, 2005 3:32 PM
> To: Gabor Grothendieck
> Cc: mingan; r-help@stat.math.ethz.ch; Huntsinger,Reid
> Subject: Re: [R] rank of a matrix
>
> Gabor Grothendieck wrote:
> > In this case, try a lower tolerance (1e-7 is the default):
> >
> >
> >>qr(hilbert(9), tol = 1e-8)$rank
> >
> > [1] 9
>
> But don't trust the results. For example, create a matrix with 4
> identical copies of hilbert(9). This still has rank 9. It's hard to
> find, though:
>
> > h9 <- hilbert(9)
> > temp <- cbind(h9, h9)
> > h9times4 <- rbind(temp, temp)
> >
> > qr(h9times4,tol=1e-7)$rank
> [1] 7
> > qr(h9times4, tol=1e-8)$rank
> [1] 10
> > qr(h9times4, tol=1e-9)$rank
> [1] 11
> > qr(h9times4, tol=1e-10)$rank
> [1] 12
>
>
> There's a tolerance that gives the right answer (1.5e-8 works for me),
> but how would I know that in a real problem where I didn't already know
> the answer?
>
> Duncan Murdoch
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-
> guide.html


R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Thu May 05 08:12:28 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:31:35 EST