[R] Solving Matrices
etb
lizzy at noradd.org
Fri Apr 16 18:09:33 CEST 2004
"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
}
More information about the R-help
mailing list