[R] Help with improving efficiency

[Bill.Venables@cmis.csiro.au]

Roger Peng suggests:

>  -----Original Message-----
> From: 	Roger Peng [mailto:rpeng at stat.ucla.edu] 
> Sent:	Saturday, August 10, 2002 2:26 PM
> To:	Wolfgang Viechtbauer
> Cc:	r-help at stat.math.ethz.ch
> Subject:	Re: [R] Help with improving efficiency
> 
> In your example, you might consider
> 
> m1 <- t(W) %*% diag(v) %*% W
> 
	[WNV]  Mathematically correct but potentially memory demanding and
slow.  What's wrong with

		m1 <- crossprod(v*W, W)

	or, if you can assume that the elements of v are non-negative:

		m1 <- crossprod(sqrt(v) * W)

	?  I have learnt to be very wary about code that uses expliticly
stored diagonal matrices.

	Bill Venables.


> -roger
> _______________________________
> UCLA Department of Statistics
> rpeng at stat.ucla.edu
> http://www.stat.ucla.edu/~rpeng
> 
> On Fri, 9 Aug 2002, Wolfgang Viechtbauer wrote:
> 
> > > The "real" answer will tell us something about R
> > > syntax, arrays, vectorization, etc. I would also love
> > > to know it.
> > >
> > > But as it happens, the function you are describing is
> > > equal to t(W)%*%W.
> > 
> > I must admit that I completely didn't notice that my code is just a
> > cumbersome way of calculating the crossproduct of the W matrix (which
> > could be accomplished even faster with the crossprod() function, as one
> > person pointed out to me off-line).
> > 
> > However, I guess I should have mentioned that the code within the loop
> > is just a subset of the operations that are being performed on the row
> > vectors.
> > 
> > So, for example, let's say that each of the k (p*p) matrices is also
> > being multiplied by a value v[i]:
> > 
> > 
> > k <- 5
> > p <- 3
> > W <- matrix(c( rep(1, k), round(rnorm(n=(p*k), mean=0, sd=1), 0)),
> >      ncol=(p+1), nrow=k)                # create a random W matrix
> > v  <- rnorm(n=k, mean=0, sd=1)		# create random v values
> > m1 <- matrix(0, nrow=p+1, ncol=p+1)
> > for (i in 1:k) {
> >    m1 <- m1 + v[i] * as.matrix(W[i,]) %*% W[i,]
> > }
> > 
> > 
> > Obviously, this is different from simlpy taking the crossproduct of W.
> > The best way to optimize the code probably depends on the types of
> > operations being performed within the loop, but I am hoping to gain a
> > general idea of how to optimize this particular code and then maybe I
> > can apply that knowledge to the actual problem I am dealing with (and
> > the "real" answer, as pointed out, might serve a nice illustration of
> > how to write efficient code in general).
> > 
> > Wolfgang
> > 
> > 
> >
> -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.
> -.-.-
> > 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
> >
> _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._.
> _._._
> > 
> 
> -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.
> -.-.-
> 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
> _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._.
> _._._
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._