very fast OLS regression?
On Wed, Mar 25, 2009 at 5:15 PM, Gavin Simpson <gavin.simpson at ucl.ac.uk> wrote:
On Wed, 2009-03-25 at 16:28 -0400, ivo welch wrote:
Dear R experts: I just tried some simple test that told me that hand computing the OLS coefficients is about 3-10 times as fast as using the built-in lm() function. ? (code included below.) ?Most of the time, I do not care, because I like the convenience, and I presume some of the time goes into saving a lot of stuff that I may or may not need. ?But when I do want to learn the properties of an estimator whose input contains a regression, I do care about speed. What is the recommended fastest way to get regression coefficients in R? ?(Is Gentlemen's weighted-least-squares algorithm implemented in a low-level C form somewhere? ?that one was always lightning fast for me.)
No one has yet mentioned Doug Bates' article in R News on this topic, which compares timings of various methods for least squares computations. Douglas Bates. Least squares calculations in R. R News, 4(1):17-20, June 2004. A Cholesky decomposition solution proved fastest in base R code, with an even faster version developed using sparse matrices and the Matrix package. You can find Doug's article here: http://cran.r-project.org/doc/Rnews/Rnews_2004-1.pdf
An updated version is available as the "Comparisons" vignette in the Matrix package.
regards,
/ivo
bybuiltin = function( y, x ) ? coef(lm( y ~ x -1 ));
byhand = function( y, x ) {
? xy<-t(x)%*%y;
? xxi<- solve(t(x)%*%x)
? b<-as.vector(xxi%*%xy)
? ## I will need these later, too:
? ## res<-y-as.vector(x%*%b)
? ## soa[i]<-b[2]
? ## sigmas[i]<-sd(res)
? b;
}
MC=500;
N=10000;
set.seed(0);
x= matrix( rnorm(N*MC), nrow=N, ncol=MC );
y= matrix( rnorm(N*MC), nrow=N, ncol=MC );
ptm = proc.time()
for (mc in 1:MC) byhand(y[,mc],x[,mc]);
cat("By hand took ", proc.time()-ptm, "\n");
ptm = proc.time()
for (mc in 1:MC) bybuiltin(y[,mc],x[,mc]);
cat("By built-in took ", proc.time()-ptm, "\n");
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
-- %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% ?Dr. Gavin Simpson ? ? ? ? ? ? [t] +44 (0)20 7679 0522 ?ECRC, UCL Geography, ? ? ? ? ?[f] +44 (0)20 7679 0565 ?Pearson Building, ? ? ? ? ? ? [e] gavin.simpsonATNOSPAMucl.ac.uk ?Gower Street, London ? ? ? ? ?[w] http://www.ucl.ac.uk/~ucfagls/ ?UK. WC1E 6BT. ? ? ? ? ? ? ? ? [w] http://www.freshwaters.org.uk %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.