Самый быстрый способ вычисления вектора ядра Гаусса в CppArmadillo? - PullRequest
Новая
       42

Самый быстрый способ вычисления вектора ядра Гаусса в CppArmadillo?

1 голос
/ 07 апреля 2019

Я пытаюсь вычислить вектор оценок ядра Гаусса как можно быстрее.У меня есть точка данных x в R ^ p и матрица X из n векторов x_i.Я хотел бы вычислить exp (- || x-x_i || ^ 2 / t) для каждого x_i и вернуть результат в виде вектора.

Я попытался реализовать это как в R, так и в RcppArmadillo через следующееcode

R CODE: 

kernel <- function(x, Data, sigma){
  if(sigma <= 0 ) stop('Gaussian kernel parameter <= 0.')

  DiffPart <- (t(t(Data) - x))^2  ## Computes the distance squared of the data and point x
  DiffPart <- rowSums(DiffPart)  # Sum of squares
  exp( - DiffPart / sigma)  #Divide by kernel parameter and evluate exponential function
}

RcppArmadillo: 

arma::Col<double> kernelCPP(arma::Row<double> x, arma::Mat<double> Data, double sigma){
  arma::Mat<double> Diff=Data.each_row()-x;
  int n = Data.n_rows;
  arma::Col<double> kern(n);
  for(int k = 0 ; k < n; k++){
    kern(k) = exp(-arma::accu(square(Diff.row(k)))/sigma);
  }
  return(kern);
}

К сожалению, мой код RcppArmadillo не намного быстрее, чем оригинальный код R.Я буду вычислять векторы ядра сотни тысяч раз в будущем коде / вычислении, и поэтому я хотел бы, чтобы этот процесс был настолько быстрым, насколько я могу.

При микробенчмаркинге я получаю следующие результаты: 

> microbenchmark(
+   kernel(x= TrainX1[1,], Data = TrainX1, sigma = 100)
+ )
Unit: milliseconds
      min       lq   mean   median
  2.223359 2.274559 2.5199 2.308052
       uq     max neval
 2.575144 4.73301   100

и 

> microbenchmark(
+   kernelCPP(x= TrainX1[1,], Data = TrainX1, sigma = 100)
+ )
Unit: milliseconds
    min       lq     mean
 1.697706 1.732053 1.826743
   median       uq      max neval
 1.775786 1.871786 2.493439   100

Чуть быстрее, но ненамного.

...