# [R] qr and Moore-Penrose

Wed Jun 30 11:12:24 CEST 1999

```yesterday I had a little shock using qr (or lm). having a matrix

X <- cbind(1,diag(3))
y <- 1:3

the qr.coef returns one NA (because X is singular). So I computed the
Moore-Penrose inverse of X (just from the definition) and I get a correct
result. Whats wrong about qr in this situation?

here is the code:

mpinv <- function(X)
{
# Moore-Penrose inverse

Eps <- 100 * .Machine\$double.eps;

# singular value decomposition

s <- svd(X);
d <- s\$d;
m <- length(d);

if (!(is.vector(d)))
return(t(s\$v%*%(1/d)%*%t(s\$u)));

# remove eigenvalues equal zero

d <- d[d > Eps];
notnull <- length(d);

if (notnull == 1)
{
inv <- 1/d;
} else {
inv <- solve(diag(d));
}

# add rows, columns of zeros if needed

if (notnull != m)
{
inv <- cbind(inv, matrix(0, nrow=notnull, ncol=(m - notnull)));
inv <- rbind(inv, matrix(0, nrow=(m-notnull), ncol=m));
}

# compute Moore-Penrose

mp <- s\$v%*%inv%*%t(s\$u);

# set very small values to zero

mp[abs(mp) < Eps] <- 0;
return(mp);
};

X <- cbind(1, diag(3));	# singular matrix
y <- 1:3
Xp <- qr(X);
b1 <- qr.coef(Xp, y);	# contains NA
b2 <- mpinv(X)%*%y	# least square fit using Moore-Penrose
X%*%b2			# == y

Torsten

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._

```