I think one needs an LU decomposition rather than QR.
However, I couldn't find anything off the shelf to do
an LU, other than learning that determinant() now
uses LU instead of QR or SVD, so the code to do it must
be in there for those that want it.
You'll probably need to divide rows of U by the first
entry if you insist on the unique reduced REF.
However, I can't see any reason not to use John Fox'
suggested R function for small and well-conditioned
matrices, such as one might have to reduce for
a homework assignment in an introductory linear algebra
course ...
Did anyone ever need a REF in any statistical application?
> -----Original Message-----
> From: John Fox [mailto:jfox at mcmaster.ca]
> Sent: 04 March 2004 15:38
> To: 'Spencer Graves'
> Cc: r-help at stat.math.ethz.ch
> Subject: RE: [R] row-echelon form (was no subject)
>
>
> Security Warning:
> If you are not sure an attachment is safe to open contact
> Andy on x234.
> There are 0 attachments with this message.
> ________________________________________________________________
>
> Dear Spencer,
>
> The R matrix from the qr decomposition isn't quite in
> row-echelon form,
> because the leading entry in each row is not 1, and the other
> entries in a
> column with a leading entry aren't all 0. Some more examples:
>
> > A # nonsingular
> [,1] [,2] [,3]
> [1,] 2 -2 0
> [2,] 1 -1 1
> [3,] 4 4 -4
> > rowEchelonForm(A)
> [,1] [,2] [,3]
> [1,] 1 0 0
> [2,] 0 1 0
> [3,] 0 0 1
> > qr.R(qr(A))
> [,1] [,2] [,3]
> [1,] -4.582576 -2.400397 3.2732684
> [2,] 0.000000 3.903600 -2.3421602
> [3,] 0.000000 0.000000 0.8944272
>
>
> > B # rank 2
> [,1] [,2] [,3] [,4]
> [1,] -2 0 -1 2
> [2,] 4 0 1 0
> [3,] 6 0 1 2
> > rowEchelonForm(B)
> [,1] [,2] [,3] [,4]
> [1,] 1 0 0 1
> [2,] 0 0 1 -4
> [3,] 0 0 0 0
> > qr.R(qr(B))
> [,1] [,2] [,3] [,4]
> [1,] 7.483315 1.6035675 0 1.069045e+00
> [2,] 0.000000 0.6546537 0 -2.618615e+00
> [3,] 0.000000 0.0000000 0 6.336077e-16
>
>
> Regards,
> John
>
> > -----Original Message-----
> > From: Spencer Graves [mailto:spencer.graves at pdf.com]
> > Sent: Thursday, March 04, 2004 9:31 AM
> > To: John Fox
> > Cc: 'Aimin Yan'; r-help at stat.math.ethz.ch
> > Subject: Re: [R] row-echelon form (was no subject)
> >
> > How about the following:
> >
> > > A <- array(1:6, dim=c(3, 2))
> > > A.qr <- qr(A)
> > > qr.R(A.qr)
> > [,1] [,2]
> > [1,] -3.741657 -8.552360
> > [2,] 0.000000 1.963961
> > >
> > I'm no expert, either, and I don't have time now to
> > research this further. Perhaps someone else can further
> > enlighten us both.
> > spencer graves
> >
> > John Fox wrote:
> >
> > >Dear Spencer,
> > >
> > >I'd be surprised if the qr decomposition as computed in R
> > weren't a lot
> > >more stable numerically, but I'm no expert. As well, I
don't
> > know how
> > >to get the row-echelon form from the qr decomposition -- though I
> > >suspect that you or someone else on the list is about to
> > enlighten me.
> > >
> > >Regards,
> > > John
> > >
> > >
> > >
> > >>-----Original Message-----
> > >>From: Spencer Graves [mailto:spencer.graves at pdf.com]
> > >>Sent: Wednesday, March 03, 2004 10:45 PM
> > >>To: John Fox
> > >>Cc: 'Aimin Yan'; r-help at stat.math.ethz.ch
> > >>Subject: Re: [R] row-echelon form (was no subject)
> > >>
> > >>How does this compare with R of the qr decomposition?
> > spencer graves
> > >>
> > >>John Fox wrote:
> > >>
> > >>
> > >>
> > >>>Dear Amin,
> > >>>
> > >>>I have a function (created just for demonstration, and
reproduced
> > >>>below) for finding the row-echelon form of a matrix.
I'm
> > >>>
> > >>>
> > >>sure that many
> > >>
> > >>
> > >>>list members could produce something that's better
> > numerically, but
> > >>>this should be OK at least for toy problems.
> > >>>
> > >>>John
> > >>>
> > >>>--------- snip -------------
> > >>>
> > >>>rowEchelonForm <- function(X, tol=.Machine$double.eps){
> > >>> if ((!is.matrix(X)) || (!is.numeric(X)))
stop("argument
> > >>>
> > >>>
> > >>must be a
> > >>
> > >>
> > >>>numeric matrix")
> > >>> Z <- X
> > >>> for (i in 1:min(dim(X))){
> > >>> if (i > 1) Z[i-1,] <- 0
> > >>> which <- which.max(abs(Z[,i])) # find maximum
pivot
> > >>>
> > >>>
> > >>in current
> > >>
> > >>
> > >>>column at or below current row
> > >>> pivot <- X[which, i]
> > >>> if (abs(pivot) <= tol) next # check for 0
pivot
> > >>> if (which > i) X[c(i,which),] <-
X[c(which,i),] #
> > >>>
> > >>>
> > >>exchange rows
> > >>
> > >>
> > >>> X[i,] <- X[i,]/pivot # pivot
> > >>> row <- X[i,]
> > >>> X <- X - outer(X[,i], row) # sweep
> > >>> X[i,] <- row # restore
current row
> > >>> }
> > >>> n <- nrow(X)
> > >>> for (i in 1:n) if (max(abs(X[i,])) <= tol)
X[c(i,n),] <-
> > >>>
> > >>>
> > >>X[c(n,i),] #
> > >>
> > >>
> > >>>0 rows to bottom
> > >>> X
> > >>> }
> > >>>
> > >>>
> > >>>
> > >>>
> > >>>
> > >>>
> > >>>>-----Original Message-----
> > >>>>From: r-help-bounces at stat.math.ethz.ch
> > >>>>[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf
Of Aimin Yan
> > >>>>Sent: Wednesday, March 03, 2004 2:42 PM
> > >>>>To: r-help at stat.math.ethz.ch
> > >>>>Subject: [R] (no subject)
> > >>>>
> > >>>>how to produce a Row Reduced Echelon Form for a
matrix in R?
> > >>>>Aimin Yan
> > >>>>
> > >>>>______________________________________________
> > >>>>R-help at stat.math.ethz.ch mailing list
> > >>>>https://www.stat.math.ethz.ch/mailman/listinfo/r-help
> > >>>>PLEASE do read the posting guide!
> > >>>>http://www.R-project.org/posting-guide.html
> > >>>>
> > >>>>
> > >>>>
> > >>>>
> > >>>______________________________________________
> > >>>R-help at stat.math.ethz.ch mailing list
> > >>>https://www.stat.math.ethz.ch/mailman/listinfo/r-help
> > >>>PLEASE do read the posting guide!
> > >>>http://www.R-project.org/posting-guide.html
> > >>>
> > >>>
> > >>>
> > >>>
> > >
> > >
> > >
> >
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://www.stat.math.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
Simon Fear
Senior Statistician
Syne qua non Ltd
Tel: +44 (0) 1379 644449
Fax: +44 (0) 1379 644445
email: Simon.Fear at synequanon.com
web: http://www.synequanon.com
Number of attachments included with this message: 0
This message (and any associated files) is confidential and\...{{dropped}}