On 03/07/12 03:58, chuck.01 wrote:

> Hi,

> I was playing around with something else and I noticed this matrix code for

> residuals in a linear model doesn't say what lm() says. Please tell me if I

> am completely misguided here.

>

> data(mtcars)

> Y <- as.matrix(mtcars[,1])

> X <- as.matrix(mtcars[,c(2:11)])

>

> # shouldnt this:

> H <- X %*% solve(t(X) %*% X) %*% t(X)

> (diag(dim(H)[1]) - H) %*% Y

>

> # be equal to this:

> residuals(lm(Y~X))

>

> # ???

> # thanks

You are forgetting about the constant term.

Try

X <- cbind(1,X)

H <- X %*% solve(t(X) %*% X) %*% t(X)

R1 <-(diag(dim(H)[1]) - H) %*% Y

R2 <-residuals(lm(Y~X))

range(R1-R2)

I get:

[1] -3.886808e-12 1.262768e-11

OMMMMMMMMMMMMMMMM!!! :-)

cheers,

Rolf Turner

______________________________________________

[hidden email] mailing list

https://stat.ethz.ch/mailman/listinfo/r-helpPLEASE do read the posting guide

http://www.R-project.org/posting-guide.htmland provide commented, minimal, self-contained, reproducible code.