Element wise array/matrix multiplication

[ace]

Suppose I have a 3 dimensional array, big_array with dimension c(5, 2, 3).

big_array <- array(
   sapply(
      rep(2, 5), mvtnorm::rmvnorm, sigma = diag(2)
   ), dim = c(5, 2, 3)
)

and a matrix, small_mat with dimension c(5, 2).

small_mat <- mvtnorm::rmvnorm(5, sigma = diag(2))

Note that the sizes of the first two dimensions coincide.

I would like to multiply them element wise so that I end up with an object with the same dimensions as big_array. What's a good way to do this?

Thanks!

[davis]

Good question! As I'm sure you've seen, you can't just do big_array * small_mat because the dimensions don't match, and R requires matching dimensions when doing array/matrix arithmetic.

One solution is to reshape small_mat to have that third dimension of "depth" 3.

library(mvtnorm)

big_array <- array(
  sapply(
    rep(2, 5), rmvnorm, sigma = diag(2)
  ), dim = c(5, 2, 3)
)

small_mat <- rmvnorm(5, sigma = diag(2))

# sigh
big_array * small_mat
#> Error in big_array * small_mat: non-conformable arrays

small_mat_reshaped <- array(small_mat, c(5, 2, 3))
small_mat_reshaped
#> , , 1
#> 
#>            [,1]       [,2]
#> [1,]  1.6497271  1.1020500
#> [2,]  0.8772834 -0.8751448
#> [3,]  0.9883818  0.5546976
#> [4,]  0.5238883 -1.1669689
#> [5,] -0.5959466  0.7387694
#> 
#> , , 2
#> 
#>            [,1]       [,2]
#> [1,]  1.6497271  1.1020500
#> [2,]  0.8772834 -0.8751448
#> [3,]  0.9883818  0.5546976
#> [4,]  0.5238883 -1.1669689
#> [5,] -0.5959466  0.7387694
#> 
#> , , 3
#> 
#>            [,1]       [,2]
#> [1,]  1.6497271  1.1020500
#> [2,]  0.8772834 -0.8751448
#> [3,]  0.9883818  0.5546976
#> [4,]  0.5238883 -1.1669689
#> [5,] -0.5959466  0.7387694

big_array * array(small_mat, c(5, 2, 3))
#> , , 1
#> 
#>             [,1]       [,2]
#> [1,] -0.27196864  0.4807641
#> [2,] -0.08161026 -1.4298487
#> [3,] -0.73607024 -0.0011445
#> [4,] -0.19076169  0.2165741
#> [5,]  0.34588880  0.9547587
#> 
#> , , 2
#> 
#>            [,1]       [,2]
#> [1,] 0.06962186 -1.6031662
#> [2,] 0.73160661  0.2134017
#> [3,] 0.31997080  0.5269535
#> [4,] 0.50023436 -1.2813516
#> [5,] 1.11332967 -1.4334219
#> 
#> , , 3
#> 
#>             [,1]       [,2]
#> [1,] -0.27196864  0.4807641
#> [2,] -0.08161026 -1.4298487
#> [3,] -0.73607024 -0.0011445
#> [4,] -0.19076169  0.2165741
#> [5,]  0.34588880  0.9547587

[davis]

As a follow up, I've been thinking about this kind of operation quite a lot as of late. Without knowing it, you're asking the question, "can R do broadcasting operations?" Broadcasting is the concept of extending dimensions when doing arithmetic. The answer is essentially no, which I find a bit frustrating. In the Python world, numpy has excellent support for this. I find broadcasting to be very intuitive, and personally think R could greatly benefit from having it.

I've begun working on a very experimental package called rray that implements this idea.

GitHub - r-lib/rray: Simple Arrays

library(mvtnorm)
library(rray)

big_array <- array(
  sapply(
    rep(2, 5), rmvnorm, sigma = diag(2)
  ), dim = c(5, 2, 3)
)

small_mat <- rmvnorm(5, sigma = diag(2))

# basically looks the same as an array
small_mat_rray <- as_rray(small_mat)
small_mat_rray
#> <vctrs_rray<double>[,2][10]>
#>      [,1]       [,2]      
#> [1,]  2.2309991 -0.9205287
#> [2,]  1.3363968 -0.2001692
#> [3,]  0.6231761  0.3299329
#> [4,] -0.4492904 -0.6025955
#> [5,]  0.7387441  1.0817863

# but broadcasting works out of the box!
big_array * small_mat_rray
#> <vctrs_rray<double>[,2,3][30]>
#> , , 1
#> 
#>       
#>        [,1]        [,2]       
#>   [1,] -4.53325710 -0.02193209
#>   [2,]  0.66145198  0.02869939
#>   [3,] -0.45921536 -0.37524633
#>   [4,] -0.51228922  0.36587511
#>   [5,] -0.88641788 -0.97357096
#> 
#> , , 2
#> 
#>       
#>        [,1]        [,2]       
#>   [1,] -2.01745243  0.35840986
#>   [2,]  1.05230941 -0.09413025
#>   [3,] -0.04026026  0.21396993
#>   [4,]  0.13486991  0.53042421
#>   [5,] -0.06176950 -1.96594265
#> 
#> , , 3
#> 
#>       
#>        [,1]        [,2]       
#>   [1,] -4.53325710 -0.02193209
#>   [2,]  0.66145198  0.02869939
#>   [3,] -0.45921536 -0.37524633
#>   [4,] -0.51228922  0.36587511
#>   [5,] -0.88641788 -0.97357096

If you want to learn more about broadcasting, I think the vignette does a good job of explaining the concept https://davisvaughan.github.io/rray/articles/broadcasting.html

[ace]

@davis has a great answer which is elegant and easy to follow - you should use that one. My answer here is not elegant or easy to follow, but I'm putting it down anyway!

big_array <- array(
  sapply(
    rep(2, 5), mvtnorm::rmvnorm, sigma = diag(2)
  ), dim = c(5, 2, 3)
)

small_mat <- mvtnorm::rmvnorm(5, sigma = diag(2))

l_array <- lapply(1:3, function(x){big_array[,,x]})
l_prod_array <- lapply(l_array, function(x, y){x * y}, y = small_mat)

prod_array <- array(sapply(l_prod_array, function(x){x}), dim = c(5, 2, 3))