Understanding the Output of Vector and Matrix Multiplication in R

76 Views Asked by At

I am currently learning about matrix operations in R, and I encountered an output that I did not anticipate when multiplying a vector with a matrix. My understanding was that multiplication would be conducted element-wise, but the output suggests otherwise. The code I ran is as follows:

num_matrix <- matrix(1:9, byrow = TRUE, nrow = 3)
new_matrix <- c(1,2,3) * num_matrix
print(new_matrix)

The output I got is:
[,1] [,2] [,3]
[1,]    1    2    3
[2,]    8   10   12
[3,]   21   24   27

This looks like the vector c(1,2,3) has been multiplied with each column of num_matrix, which contradicts my initial expectation of a row-wise operation. Could someone explain why the multiplication behaves like this and not row-wise as I expected?

Additionally, if I want to multiply a vector with each row of a matrix, rather than with each column, how can I achieve that in R?

Thank you for any assistance!

2

There are 2 best solutions below

1
PBulls On

Matrices in R are vectors in column-major order. You can transpose your matrix to have elements recycled along columns that were previously rows, and then transpose back:

(want_matrix <- t(t(num_matrix) * 1:3))
#>      [,1] [,2] [,3]
#> [1,]    1    4    9
#> [2,]    4   10   18
#> [3,]    7   16   27

For kicks, here's a few other methods and their benchmarks. Note that treating this as a matrix multiplication will be fastest when the input is small, but doesn't scale well -- transposing is more consistent.

Edit: the method of @ThomasIsCoding (f5 below) turns out to be faster throughout!

f1 <- \(m, v) t(t(m) * v)
f2 <- \(m, v) t(apply(m, 1, \(r) r*v))
f3 <- \(m, v) sweep(m, 2, v, `*`)
f4 <- \(m, v) m %*% diag(v)
f5 <- \(m, v) m * v[col(m)]

bench <- function(n, ...) {
  m <- matrix(rnorm(n*n), ncol=n)
  v <- runif(n)

  microbenchmark::microbenchmark(
    transpose=f1(m,v),
    apply=f2(m,v),
    sweep=f3(m,v),
    `%*%`=f4(m,v),
    col=f5(m,v),
    ...,
    check = "equivalent"
  )
}

bench(10)
#> Unit: microseconds
#>       expr   min     lq    mean median     uq   max
#>  transpose  30.4  33.40  34.919   34.1  35.00  53.1
#>      apply 121.5 124.60 128.871  125.6 128.60 188.1
#>      sweep  94.1  97.90 107.103   99.5 102.80 275.0
#>        %*%  13.1  14.70  15.272   15.1  15.45  29.9
#>        col  10.9  12.05  12.548   12.5  12.90  21.0

bench(100)
#> Unit: microseconds
#>       expr   min     lq    mean median     uq   max
#>  transpose  92.0  95.50 101.730  96.90 102.45 157.6
#>      apply 392.3 398.05 418.783 403.60 427.35 566.8
#>      sweep 154.4 159.40 171.652 161.20 168.00 390.7
#>        %*% 495.4 497.20 514.902 498.00 524.75 907.0
#>        col  60.7  62.15  65.206  62.95  65.90 117.3
0
ThomasIsCoding On

You can try col+ * like below and benchmark

microbenchmark(
    f1 = num_matrix * (1:3)[col(num_matrix)],
    f2 = num_matrix %*% diag(1:3),
    unit = "relative",
    check = "equivalent"
)

and you will see

Unit: relative
 expr      min       lq    mean median       uq       max neval
   f1 1.000000 1.000000 1.00000 1.0000 1.000000 1.0000000   100
   f2 1.553829 1.498501 1.30164 1.4995 1.454133 0.5962642   100