У меня есть матрица mat и вектор v. Я хотел бы умножить первый столбец матрицы mat на первый элемент вектора v и умножить второй столбец матрицы mat на второй элемент вектора v. Я могу сделать это, как показано. Как я могу сделать это быстрее в R, так как мы получаем большую матрицу?
mat = matrix(rnorm(1500000), ncol= 100)
v= rnorm(100)
> system.time( mat %*% diag(v))
user system elapsed
0.02 0.00 0.02
sweep кажется, работает немного быстрее на моей машине
sweep(mat, 2, v, FUN = "*")
Некоторые ориентиры:
> microbenchmark(mat %*% diag(v),sweep(mat, 2, v, FUN = "*"))
Unit: milliseconds
expr min lq median uq max neval
%*% 214.66700 226.95551 231.2366 255.78493 349.1911 100
sweep 42.42987 44.72254 62.9990 70.87403 127.2869 100
Утилизация может ускоряться, но вы перерабатываете в столбцах, а не поперек, поэтому просто переносите и переносите назад.
t( t(mat) * v )
Это должно быть быстрее, чем sweep или %*%.
microbenchmark(mat %*% diag(v),sweep(mat, 2, v, FUN = "*"), t(t(mat)*v))
Unit: milliseconds
expr min lq median uq max neval
%*% 150.47301 152.16306 153.17379 161.75416 281.3315 100
sweep 35.94029 42.67210 45.53666 48.07468 168.3728 100
t(t(mat) * v) 16.50813 23.41549 26.31602 29.44008 160.1651 100
Немного поздно в игре, но кто-то сказал быстрее всего? Это может быть другим хорошим использованием Rcpp. Эта функция (называемая mmult) будет по умолчанию умножать каждый столбец матрицы на каждый последующий элемент вектора, но имеет возможность сделать это по столбцу, установив byrow = FALSE. Он также проверяет, что m и v имеют соответствующий размер с опцией byrow. Во всяком случае, он быстрый (примерно в 10-12 раз быстрее, чем лучший ответ на родной R)...
@chris предоставил этот отличный ответ к другому вопросу, который я задал, пытаясь заставить это работать с RcppArmadillo. Однако кажется, что функция Rcpp -одно, которую я размещал здесь, по-прежнему примерно в 8 раз быстрее, чем примерно, и примерно в 70 раз быстрее, чем метод OP. Нажмите ссылку для кода для функции @chris '- это красиво просто.
Я поставлю бенчмаркинг наверху.
require( microbenchmark )
m <- microbenchmark( mat %*% diag(v) , mmult( mat , v ) , sweep(mat, 2, v, FUN = "*") , chris( mat , v ) , t( t(mat) * v ) , times = 100L )
print( m , "relative" , order = "median" , digits = 3 )
Unit: relative
expr min lq median uq max neval
mmult(mat, v) 1.00 1.00 1.00 1.00 1.00 100
chris(mat, v) 10.74 9.31 8.15 7.27 10.44 100
t(t(mat) * v) 9.65 8.75 8.30 15.33 9.52 100
sweep(mat, 2, v, FUN = "*") 20.51 18.35 22.18 21.39 16.94 100
mat %*% diag(v) 80.44 70.11 73.12 70.68 54.96 100
Просмотрите, как работает mmult и возвращает тот же результат, что и OP...
require( Rcpp )
# Source code for our function
func <- 'NumericMatrix mmult( NumericMatrix m , NumericVector v , bool byrow = true ){
if( byrow );
if( ! m.nrow() == v.size() ) stop("Non-conformable arrays") ;
if( ! byrow );
if( ! m.ncol() == v.size() ) stop("Non-conformable arrays") ;
NumericMatrix out(m) ;
if( byrow ){
for (int j = 0; j < m.ncol(); j++) {
for (int i = 0; i < m.nrow(); i++) {
out(i,j) = m(i,j) * v[j];
}
}
}
if( ! byrow ){
for (int i = 0; i < m.nrow(); i++) {
for (int j = 0; j < m.ncol(); j++) {
out(i,j) = m(i,j) * v[i];
}
}
}
return out ;
}'
# Make it available
cppFunction( func )
# Use it
res1 <- mmult( m , v )
# OP function
res2 <- mat %*% diag(v)
# Same result?
identical( res1 , res2 ) # Yes!!
[1] TRUE