оптимизировать матричное умножение для l oop RcppArmadillo - PullRequest
Новая
       78

оптимизировать матричное умножение для l oop RcppArmadillo

0 голосов
/ 20 июня 2020

Цель состоит в том, чтобы реализовать быструю версию ортогональной проективной неотрицательной матричной факторизации (opnmf) в R. Я перевожу имеющийся код MATLAB здесь .

Я реализовал версия vanilla R, но она намного медленнее (примерно в 5,5 раза медленнее), чем реализация matlab на моих данных (~ 225000 x 150) для 20-факторного решения.

Итак, я подумал, что использование C ++ может ускорить работу, но его скорость похож на R. Я думаю, что это можно оптимизировать, но не уверен, как это сделать, поскольку я новичок ie в c ++. Здесь - это поток, в котором обсуждается аналогичная проблема.

Вот моя реализация RcppArmadillo. 

// [[Rcpp::export]]
Rcpp::List arma_opnmf(const arma::mat & X, const arma::mat & W0, double tol=0.00001, int maxiter=10000, double eps=1e-16) {
  arma::mat W = W0;
  arma::mat Wold = W;
  arma::mat XXW = X * (X.t()*W);
  double diffW = 9999999999.9;
  
  Rcout << "The value of maxiter : " << maxiter << "\n";
  Rcout << "The value of tol : " << tol << "\n";
  
  int i;
  for (i = 0; i < maxiter; i++) {
    XXW = X * (X.t()*W);
    W = W % XXW / (W  * (W.t() * XXW));
    //W = W % (X*(X.t()*W)) / (W*((W.t()*X)*(X.t()*W)));
    
    arma::uvec idx = find(W < eps);
    W.elem(idx).fill(eps);
    W = W / norm(W,2);
    diffW = norm(Wold-W, "fro") / norm(Wold, "fro");
    if(diffW < tol) {
      break;
    } else {
      Wold = W;
    }
    
    if(i % 10 == 0) {
      Rcpp::checkUserInterrupt();
    }
    
  }
  return Rcpp::List::create(Rcpp::Named("W")=W,
                            Rcpp::Named("iter")=i,
                            Rcpp::Named("diffW")=diffW);
}

Эта предлагаемая проблема подтверждает, что matlab является довольно быстро, так что при использовании R / c ++ надежды нет?

Тесты проводились на Windows 10 и Ubuntu 16 с R версии 4.0.0.

EDIT

После интересных комментариев в ответе ниже. Размещаю дополнительные подробности. Я провел тесты на машине Windows 10 с R 3.5.3 (как это и предлагает Microsoft), и сравнение показывает, что RcppArmadillo с Microsoft R является самым быстрым.

R 

   user  system elapsed 
 213.76    7.36  221.42

R с RcppArmadillo 

   user  system elapsed 
 179.88    3.44  183.43

Microsoft Open R 

   user  system elapsed 
 167.33    9.96   45.94

Microsoft Open с RcppArmadillo 

    user  system elapsed 
  85.47    4.66   23.56

1 Ответ

3 голосов
/ 20 июня 2020

Знаете ли вы, что этот код «в конечном итоге» выполняется парой библиотек, называемых LAPACK и BLAS?

Знаете ли вы, что Matlab поставляется с высоко оптимизированной библиотекой? Знаете ли вы, что во всех системах, на которых работает R, вы можете изменить , какой LAPACK / BLAS используется.

Разница имеет значение очень . Буквально сегодня утром друг опубликовал этот твит , противопоставляющий тот же код R , работающий на том же Windows компьютере , но в двух разных средах R. в шесть раз быстрее один «просто» использует параллельную реализацию LAPACK / BLAS.

Здесь вы даже не сказали нам, в какой операционной системе вы работаете. Вы можете получить OpenBLAS (который использует параллелизм) для всех ОС, на которых работает R. Вы даже можете довольно легко получить Intel MKL (который IIR C также использует Matlab) в некоторых ОС. Для Ubuntu / Debian я опубликовал скрипт на GitHub , который делает это за один шаг.

Наконец, много лет go я «унаследовал» быструю программу, работающую в Matlab на ( then-large-i sh) Windows компьютер. Я переписал часть Matlab (осторожно и медленно, это усилие) на C ++ с использованием RcppArmadillo, что привело к некоторым улучшениям - и поскольку мы могли запускать этот (теперь открытый исходный код) код параллельно из R на том же компьютере, еще несколько факторов. Вместе это было на несколько порядков превратить дневную симуляцию во что-то, что длилось несколько минут. Так что «да, вы можете».

Изменить: Поскольку у вас есть доступ к Ubuntu, вы можете переключиться с basi c LAPACK / BLAS на OpenBLAS с помощью одной команды, хотя я больше не знаком с Ubuntu 16.04 (поскольку я сам запускаю 20.04).

Edit 2: Взяв сравнение из твита Йозефа , Docker * 1033 Контейнер * r-base Я также поддерживаю (как часть Rocker Project ) OpenBLAS. [1] Итак, как только мы добавим его, например, через apt-get install libopenblas-dev, время простого повторяющегося матричного перекрестного произведения переместится с 

root@0eb44b1fcc06:/# Rscript -e 'v <- matrix(1:1e6,1e3); system.time(replicate(10, crossprod(v,v)))'
   user  system elapsed 
  9.289   0.084   9.373 
root@0eb44b1fcc06:/#

на 

root@67bd334f53d4:/# Rscript -e 'v <- matrix(1:1e6,1e3); system.time(replicate(10, crossprod(v,v)))'
   user  system elapsed 
  2.259   2.370   0.447 
root@67bd334f53d4:/#

, что является существенным.

...