[R] Solving Matrices
John Fox
jfox at mcmaster.ca
Fri Apr 16 19:59:59 CEST 2004
Dear Elizabeth,
If this is the rowEchelonForm function that I posted to r-help some
time ago, as I mentioned in a subsequent post, the tolerance was set
too low; if you set the tol argument (better, the default) to a larger
number than the machine precision, then the function should work more
reliably.
I hope this helps,
John
On 16 Apr 2004 11:09:33 -0500
etb <lizzy at noradd.org> wrote:
> "Richard" writes:
>
> [ ... ]
> >
> > So here are three different functions all saying much the same
> thing
> > about x: it is numerically close to a matrix which does NOT span
> the
> > whole of R**4.
>
> Theorem 4 (in my text):
>
> Let A be an m x n matrix. Then the following statements are
> logically equivalent. That is, for a particular A, either they
> are all true statements or they are all false:
>
> a) For each b in R**m, the equation Ax = b has a solution.
> b) The columns of A span R**m.
> c) A has a pivot position in every row.
>
> The back of the book is in agreement with you but I can manipulate
> the
> matrix to produce (c) above which should mean that (b) is correct:
>
> > x <- matrix(data=c(7,-5,6,-7,2,-3,10,9,-5,4,-2,2,8,-9,7,15),nrow=4)
> > rowEchelonForm(x)
> [,1] [,2] [,3] [,4]
> [1,] 1 0 0 0
> [2,] 0 1 0 0
> [3,] 0 0 1 0
> [4,] 0 0 0 1
>
> Likewise, I can produce the above using the swap(), scale(), gauss()
> and bgauss() functions that the author of the text originally wrote
> for maple, Mathematica, etc. but from all I can tell, the row echelon
> form of matrix x meets the criteria of theorem 4 above, yet the
> columns
> of x do not span R**4.
>
> Using the gauss()[1] function:
>
> [,1] [,2] [,3] [,4]
> [1,] 7 2.000000 -5.0000000 8.000000e+00
> [2,] 0 -1.571429 0.4285714 -3.285714e+00
> [3,] 0 0.000000 4.5454545 -1.718182e+01
> [4,] 0 0.000000 0.0000000 -5.035972e-15
>
> Again, this has a pivot column in each row although [4,4] seems a bit
> suspicious.
>
> Elizabeth
>
> [1]:
>
> gauss <- function(A, r, v = c()) {
> m <- dim(A)[1]
> n <- dim(A)[2]
> if(length(r) > 1 | r < 1 | r > m)
> stop("The second entry in gauss(...) must be the pivot row
> index.")
> if(length(v) == 0) {
> if(r+1 > m)
> v <- m
> else
> v <- (r+1):m
> }
> for(j in v) {
> if(j == r)
> stop("Pivot ow cannot change itself.")
> }
>
> col <- 1
> while(A[r, col] == 0 && col < n)
> col <- col + 1
> if(col == n && A[r, col] == 0)
> stop("Row r has only zeroes and cannot be a pivot row.")
> for(j in v) {
> multiplier <- -A[j, col]/A[r, col]
> A[j,] <- A[j,] + multiplier * A[r,]
> }
> A
> }
>
> ______________________________________________
> 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
--------------------------------
John Fox
Department of Sociology
McMaster University
Hamilton, Ontario, Canada
http://socserv.mcmaster.ca/jfox/
More information about the R-help
mailing list