multiplying multidimensional arrays (was: Re: [R] Manipulation involving arrays)

Sun, Jul 16, 2006 7:29 PM #

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 mulitplied
along 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 at gmail.com> wrote:

Ravi Varadhan

Mon, Jul 17, 2006 6:18 AM #

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 at 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 at 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 at 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]
 }