From: RAVI VARADHAN <rvaradhan_at_jhmi.edu>

Date: Mon 17 Jul 2006 - 13:18:45 GMT

R-devel@r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-devel Received on Mon Jul 17 23:22:34 2006

Date: Mon 17 Jul 2006 - 13:18:45 GMT

Dear Gabor,

Thank you very much for your solution. It speeded calculations by a factor of two.

Now to your recommendation about making this solution a basic operation. I completely agree with your suggestion. That is exactly what I would have hoped for. In fact, my first try was to do:

Tbar <- tarray %*% t(wt)

Since tarray is (3,3,mcsamp) and wt is (n,mcsamp), I figured that the result would be a (3,3,n) array that would sum over mcsamp, as in the case of 2-dim array multiplication. But obviously that didn't work as "non-conformable". Your simple trick of forcing the 3-dim array to a 2-dim matrix using matrix(tarray, 3*3) is quite clever.

Thanks again for your help.

Best,

Ravi.

- Original Message ----- From: Gabor Grothendieck <ggrothendieck@gmail.com> Date: Sunday, July 16, 2006 10:29 pm Subject: multiplying multidimensional arrays (was: Re: [R] Manipulation involving arrays)

> I am moving this to r-devel.

*>
**> The problem and solution below posted on r-help could have been
**> a bit slicker if %*% worked with multidimensional arrays multiplying
**> them so that if the first arg is a multidimensional array it is
**> mulitpliedalong the last dimension (and first dimension for the
**> second arg).
**> Then one could have written:
**>
**> Tbar <- tarray %*% t(wt) / rep(wti, each = 9)
**>
**> which is a bit nicer than what had to be done, see below, given
**> that %*% only
**> works with matrices.
**>
**> I suggest that %*% be so extended to multidimensional arrays. Note
**> that this is upwardly compatible and all existing cases would continue
**> to work unchanged.
**>
**> On 7/16/06, Gabor Grothendieck <ggrothendieck@gmail.com> wrote:
**> > The double loop is the same as:
**> >
**> > Tbar[] <- matrix(tarray, 9) %*% t(wt) / rep(wti, each = 9)
**> >
**> >
**> > On 7/16/06, RAVI VARADHAN <rvaradhan@jhmi.edu> wrote:
**> > > Hi,
**> > >
**> > > I have the following piece of code that is part of a larger
**> function. This piece is the most time consuming part of the
**> function, and I would like to make this a bit more efficient.
**> Could anyone suggest a way to do this faster?
**> > >
**> > > In particular, I would like to replace the nested "for" loop
**> with a faster construct. I tried things like "kronecker" and
**> "outer" combined with apply, but couldn't get it to work.
**> > >
**> > >
**> > > Here is a sample code:
**> > >
**> > > ##########################
**> > > n <- 120
**> > > sigerr <- 5
**> > > covmat <- diag(c(8,6,3.5))
**> > > mu <- c(105,12,10)
**> > > mcsamp <- 10000
**> > >
**> > > Tbar <- array(0, dim=c(3,3,n))
**> > >
**> > > # theta is a mcsamp x 3 matrix
**> > > theta <- mvrnorm(mcsamp, mu = mu, Sigma = covmat)
**> > >
**> > > wt <- matrix(runif(n*mcsamp),n,mcsamp)
**> > > wti <- apply(wt,1,sum)
**> > >
**> > > tarray <-
**> array(apply(theta,1,function(x)outer(x,x)),dim=c(3,3,mcsamp))> >
**> > > for (i in 1:n) {
**> > > for (k in 1:mcsamp) {
**> > > Tbar[,,i] <- Tbar[,,i] + wt[i,k] * tarray[,,k]
**> > > }
**> > > Tbar[,,i] <- Tbar[,,i] / wti[i]
**> > > }
**> > >
**> >
*

>

R-devel@r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-devel Received on Mon Jul 17 23:22:34 2006

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.1.8, at Mon 17 Jul 2006 - 22:31:08 GMT.

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help.
Please read the posting
guide before posting to the list.
*