Самый быстрый способ умножения матричных столбцов с векторными элементами в R
У меня есть матрица m
и вектор v
. Я хотел бы умножить первый столбец матрицы m
на первый элемент вектора v
и умножить второй столбец матрицы m
на второй элемент вектора v
и так далее. Я могу сделать это со следующим кодом, но я ищу способ, который не требует двух транспозиционных вызовов. Как я могу сделать это быстрее в R?
m <- matrix(rnorm(120000), ncol=6)
v <- c(1.5, 3.5, 4.5, 5.5, 6.5, 7.5)
system.time(t(t(m) * v))
# user system elapsed
# 0.02 0.00 0.02
Ответы
Ответ 1
Используйте некоторую линейную алгебру и выполняйте матричное умножение, которое довольно быстро выполняется в R
.
eg
m %*% diag(v)
некоторый бенчмаркинг
m = matrix(rnorm(1200000), ncol=6)
v=c(1.5, 3.5, 4.5, 5.5, 6.5, 7.5)
library(microbenchmark)
microbenchmark(m %*% diag(v), t(t(m) * v))
## Unit: milliseconds
## expr min lq median uq max neval
## m %*% diag(v) 16.57174 16.78104 16.86427 23.13121 109.9006 100
## t(t(m) * v) 26.21470 26.59049 32.40829 35.38097 122.9351 100
Ответ 2
Если у вас больше столбцов, ваше решение t (t (m) * v) превосходит решение по умножению матриц с большим отрывом. Тем не менее, есть более быстрое решение, но оно требует больших затрат на использование памяти. Вы создаете матрицу размером m с помощью rep() и умножаете поэлементно. Вот сравнение, модифицирующее пример mnel:
m = matrix(rnorm(1200000), ncol=600)
v = rep(c(1.5, 3.5, 4.5, 5.5, 6.5, 7.5), length = ncol(m))
library(microbenchmark)
microbenchmark(t(t(m) * v),
m %*% diag(v),
m * rep(v, rep.int(nrow(m),length(v))),
m * rep(v, rep(nrow(m),length(v))),
m * rep(v, each = nrow(m)))
# Unit: milliseconds
# expr min lq mean median uq max neval
# t(t(m) * v) 17.682257 18.807218 20.574513 19.239350 19.818331 62.63947 100
# m %*% diag(v) 415.573110 417.835574 421.226179 419.061019 420.601778 465.43276 100
# m * rep(v, rep.int(nrow(m), ncol(m))) 2.597411 2.794915 5.947318 3.276216 3.873842 48.95579 100
# m * rep(v, rep(nrow(m), ncol(m))) 2.601701 2.785839 3.707153 2.918994 3.855361 47.48697 100
# m * rep(v, each = nrow(m)) 21.766636 21.901935 23.791504 22.351227 23.049006 66.68491 100
Как видите, использование "each" в rep() жертвует скоростью для ясности. Разница между rep.int и rep кажется незначительной, обе реализации меняются местами при повторных запусках microbenchmark(). Имейте в виду, что ncol (m) == length (v).
![autoplot]()
Ответ 3
Как отмечает @Arun, я не знаю, что вы будете бить свое решение с точки зрения эффективности времени. Что касается понятности кода, есть и другие варианты:
Один вариант:
> mapply("*",as.data.frame(m),v)
V1 V2 V3
[1,] 0.0 0.0 0.0
[2,] 1.5 0.0 0.0
[3,] 1.5 3.5 0.0
[4,] 1.5 3.5 4.5
И еще:
sapply(1:ncol(m),function(x) m[,x] * v[x] )
Ответ 4
Ради полноты я добавил sweep
в тест. Несмотря на несколько вводящие в заблуждение имена атрибутов, я думаю, что он может быть более читабельным, чем другие альтернативы, а также довольно быстрым:
n = 1000
M = matrix(rnorm(2 * n * n), nrow = n)
v = rnorm(2 * n)
microbenchmark::microbenchmark(
M * rep(v, rep.int(nrow(M), length(v))),
sweep(M, MARGIN = 2, STATS = v, FUN = '*'),
t(t(M) * v),
M * rep(v, each = nrow(M)),
M %*% diag(v)
)
Unit: milliseconds
expr min lq mean
M * rep(v, rep.int(nrow(M), length(v))) 5.259957 5.535376 9.994405
sweep(M, MARGIN = 2, STATS = v, FUN = '*') 16.083039 17.260790 22.724433
t(t(M) * v) 19.547392 20.748929 29.868819
M * rep(v, each = nrow(M)) 34.803229 37.088510 41.518962
M %*% diag(v) 1827.301864 1876.806506 2004.140725
median uq max neval
6.158703 7.606777 66.21271 100
20.479928 23.830074 85.24550 100
24.722213 29.222172 92.25538 100
39.920664 42.659752 106.70252 100
1986.152972 2096.172601 2432.88704 100
Ответ 5
Как сделано bluegrue, достаточно простого представления, чтобы выполнить умножение по элементам.
Количество умножений и суммирования уменьшается на широкий диапазон, как если бы выполнялось простое матричное умножение с diag()
, где для этого случая можно было бы избежать большого количества нулевых умножений.
m = matrix(rnorm(1200000), ncol=6)
v=c(1.5, 3.5, 4.5, 5.5, 6.5, 7.5)
v2 <- rep(v,each=dim(m)[1])
library(microbenchmark)
microbenchmark(m %*% diag(v), t(t(m) * v), m*v2)
Unit: milliseconds
expr min lq mean median uq max neval cld
m %*% diag(v) 11.269890 13.073995 16.424366 16.470435 17.700803 95.78635 100 b
t(t(m) * v) 9.794000 11.226271 14.018568 12.995839 15.010730 88.90111 100 b
m * v2 2.322188 2.559024 3.777874 3.011185 3.410848 67.26368 100 a
Ответ 6
Разреженные матрицы
Если вы работаете с разреженными матрицами, вероятно, что поэлементные методы (т.е. расширение v
до плотной матрицы) могут превысить ограничения памяти, поэтому я собираюсь исключить эти предлагаемые методы из этого анализа.
Пакет Matrix предоставляет метод Diagonal
, который значительно улучшает метод диагонального умножения.
Ориентиры
library(Matrix)
library(microbenchmark)
N_ROW = 5000
N_COL = 10000
M <- rsparsematrix(N_ROW, N_COL, density=0.1)
v <- rnorm(N_COL)
microbenchmark(
M %*% Diagonal(length(v), v),
t(t(M) * v),
sweep(M, MARGIN=2, STATS=v, FUN='*'))
# Unit: milliseconds
# expr min lq mean median uq max neval
# M %*% Diagonal(length(v), v) 36.46755 39.03379 47.75535 40.41116 43.30241 248.1036 100
# t(t(M) * v) 207.70560 223.35126 269.08112 250.25379 284.49382 511.0403 100
# sweep(M, MARGIN = 2, STATS = v, FUN = "*") 1682.45263 1869.87220 1941.54691 1924.80218 1999.62484 3104.8305 100