R is performing OLS alright, the problem is in the example you provide. Here is a demonstration, building on what @DWin has commented.
set.seed(1234)
x1 <- rnorm(5,mean=3)
x2 <- rnorm(5,mean=1,sd=5)
x3 <- rnorm(5,mean=7,sd=1)
x4 <- rnorm(5,mean=1,sd=2)
X<-as.matrix(cbind(x1,x2,x3,x4))
Y<-as.matrix(cbind(y))
(beta.hat<-solve(t(X)%*%X)%*%t(X)%*%Y)
lm(Y~X+0)
As you can see, the coefficients are exactly the same, and your code is repeated except for the more suitable data.
P.S. If this is classified as an off-topic question, feel free to delete the answer as well. My intention was just to illustrate the issue with some code, which doesn't fit in a commentary.