You might use lsfit instead and just do the whole Y matrix at once. That saves all the recalculation of things involving only X.
Reid Huntsinger -----Original Message----- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Eduardo Leoni Sent: Thursday, March 03, 2005 5:16 PM To: [email protected] Subject: [R] regression on a matrix Hi - I am doing a monte carlo experiment that requires to do a linear regression of a matrix of vectors of dependent variables on a fixed set of covariates (one regression per vector). I am wondering if anyone has any idea of how to speed up the computations in R. The code follows: #regression function #Linear regression code qreg <- function(y,x) { X=cbind(1,x) m<-lm.fit(y=y,x=X) p<-m$rank r <- m$residuals n <- length(r) rss <- sum(r^2) resvar <- rss/(n - p) Qr <- m$qr p1 <- 1:p R <- chol2inv(Qr$qr[p1, p1, drop = FALSE]) se <- sqrt(diag(R) * resvar) b <- m$coefficients return(c(b[2],se[2])) } #simulate a <- c(1,.63,.63,1) a <- matrix(a,2,2) c <- chol(a) C <- 0.7 theta <- 0.8 sims <- 1000 n<-20 u <- rnorm(n,0,sqrt(1-C)) w <- rgamma(n,C/theta,1/theta) e <- rnorm(n,0,sqrt(w)) x1 <- rnorm(n) x <- x1*c[2,2]+c[1,2]*w v <- e+u y <- 1+x+v w <- rgamma(n,C/theta,1/theta) #create matrix of dep variable newdep <- matrix(rnorm(length(y)*sims,y,sqrt(w)),c(length(y),sims)) monte <- apply(newdep,2,qreg,x=x) ______________________________________________ [email protected] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html ______________________________________________ [email protected] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
