Блог

Главная - Вопросы по программированию - Как я могу воспроизвести функциональность R с умножением матрицы на вектор по элементам в Rcpp или Armadillo?

Как я могу воспроизвести функциональность R с умножением матрицы на вектор по элементам в Rcpp или Armadillo?

#c #r #matrix #rcpp #rcpparmadillo

Вопрос:

В R умножение матрицы на вектор по умолчанию выполняется по элементам и работает следующим образом:

 A <- matrix(c(1,2,3,4,5,6), nrow=2)
b <- c(1,3)
A * b
     [,1] [,2] [,3]
[1,]    1    3    5
[2,]    6   12   18

По существу, каждый элемент матрицы умножается на элемент векторного циклического перебора. Я хочу воссоздать эту функциональность в C , используя Rcpp или его расширения, такие как RcppEigen или RcppArmadillo.

С Броненосцем я попробовал следующее:

 arma::mat X ;
arma::vec beta ;

beta.resize ( 2 ) ;

beta (0) = 1.0 ;
beta (1) = 3.0 ;

X.resize ( 3, 2 ) ;

X (0,0) = 1.0 ;
X (0,1) = 2.0 ;
X (1,0) = 3.0 ;
X (1,1) = 4.0 ;
X (2,0) = 5.0 ;
X (2,1) = 6.0 ;

Rcout << X * beta << std::endl ;

Что приводит к вектору: [7, 15, 23], как если бы это было матричное умножение. Замена * на a % для поэлементного умножения приводит к ошибке: «ошибка: поэлементное умножение: несовместимые размеры матрицы: 3×2 и 2×1»

Есть ли какая-либо встроенная функция для такого R-подобного поведения, или мне нужно будет создать свою собственную функцию (которая у меня есть, но я бы предпочел «официальный» способ сделать это)

 NumericMatrix matrixTimesVector(NumericMatrix mat, NumericVector vec){
    NumericMatrix res(mat.rows(), mat.cols());
    int index = 0;
    for (int i = 0; i < mat.cols();   i){
        for (int j = 0; j < mat.rows();   j){
            res( j , i ) = mat( j , i ) * vec[index   % vec.size()];
        }
    }
    return res;
}

Ответ №1:

Вот несколько вариантов. Я полагаю, что третий ближе всего к тому, что вы хотите, — задокументирован в http://arma.sourceforge.net/docs.html#each_colrow

 #include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]

// looping through each column and element wise multiplication
// [[Rcpp::export]]
arma::mat matTimesVec(arma::mat mat, arma::vec v) {
  for(int i; i < mat.n_cols; i  ){
    mat.col(i)  %=  v;
  }
  return mat;
}

// form a diagonal matrix with the vector and then use matrix multiplication
// [[Rcpp::export]]
arma::mat matTimesVec2(arma::mat mat, arma::vec v) {
  return arma::diagmat(v) * mat;
}

// use the functionality described at http://arma.sourceforge.net/docs.html#each_colrow 
// to "Apply a vector operation to each column or row of a matrix "
// [[Rcpp::export]]
arma::mat matTimesVec3(arma::mat mat, arma::vec v) {
  mat.each_col() %= v;
  return mat; 
}

/*** R
A <- matrix(c(1,2,3,4,5,6), nrow=2)
b <- c(1,3)
A * b
matTimesVec(A, b)
matTimesVec2(A, b)
matTimesVec3(A, b)
*/
Метки: