Skip to content

(no subject)

6 messages · Aimin Yan, John Fox, Spencer Graves

Original
Aimin Yan
# 
how to produce a  Row Reduced Echelon Form for a matrix in R?
Aimin Yan
John Fox
# 
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
    }
Spencer Graves
# 
How does this compare with R of the qr decomposition?  spencer graves
John Fox wrote:

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

      
      
    

  
  


  


      

    

    

      
      

        


  

        

      John Fox
    

    

      
      #
    

  


  

      
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

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

      
      
    

            

      
      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Spencer Graves
    

    

      
      #
    

  


  

      
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:

        

            

      
      
      
    

        

      
      
    

        

      
      
    

            

      
      
      
    

  
  


  


      

    

    

      
      

        


  

        

      John Fox
    

    

      
      #
    

  


  

      
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:

            

      
      
      
    

        
[,1] [,2] [,3]
[1,]    2   -2    0
[2,]    1   -1    1
[3,]    4    4   -4

            

      
      
      
    

        
[,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1

            

      
      
      
    

        
[,1]      [,2]       [,3]
[1,] -4.582576 -2.400397  3.2732684
[2,]  0.000000  3.903600 -2.3421602
[3,]  0.000000  0.000000  0.8944272

            

      
      
      
    

        
[,1] [,2] [,3] [,4]
[1,]   -2    0   -1    2
[2,]    4    0    1    0
[3,]    6    0    1    2

            

      
      
      
    

        
[,1] [,2] [,3] [,4]
[1,]    1    0    0    1
[2,]    0    0    1   -4
[3,]    0    0    0    0

            

      
      
      
    

        
[,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