I don't have enough rep to comment, but the function given above by soldier.moth in the accepted answer [edit 2018: no longer the accepted answer] is buggy - it doesn't handle matrices where the RREF solution has zeroes on its main diagonal. Try e.g.
m
并注意输出不在 RREF 中。
我想我已经可以正常工作了,但您可能想自己检查输出:
rref <- function(A, tol=sqrt(.Machine$double.eps),verbose=FALSE,
fractions=FALSE){
## A: coefficient matrix
## tol: tolerance for checking for 0 pivot
## verbose: if TRUE, print intermediate steps
## fractions: try to express nonintegers as rational numbers
## Written by John Fox
# Modified by Geoffrey Brent 2014-12-17 to fix a bug
if (fractions) {
mass <- require(MASS)
if (!mass) stop("fractions=TRUE needs MASS package")
}
if ((!is.matrix(A)) || (!is.numeric(A)))
stop("argument must be a numeric matrix")
n <- nrow(A)
m <- ncol(A)
x.position<-1
y.position<-1
# change loop:
while((x.position<=m) & (y.position<=n)){
col <- A[,x.position]
col[1:n < y.position] <- 0
# find maximum pivot in current column at or below current row
which <- which.max(abs(col))
pivot <- col[which]
if (abs(pivot) <= tol) x.position<-x.position+1 # check for 0 pivot
else{
if (which > y.position) { A[c(y.position,which),]<-A[c(which,y.position),] } # exchange rows
A[y.position,]<-A[y.position,]/pivot # pivot
row <-A[y.position,]
A <- A - outer(A[,x.position],row) # sweep
A[y.position,]<-row # restore current row
if (verbose)
if (fractions) print(fractions(A))
else print(round(A,round(abs(log(tol,10)))))
x.position<-x.position+1
y.position<-y.position+1
}
}
for (i in 1:n)
if (max(abs(A[i,1:m])) <= tol)
A[c(i,n),] <- A[c(n,i),] # 0 rows to bottom
if (fractions) fractions (A)
else round(A, round(abs(log(tol,10))))
}