Semejanza cosina rápida de dos matrices en Rcpp con Armadillo

Estoy tratando de pasar una función R muy rápida para calcular la similitud cosina en Rcpp con Armadillo y operaciones de matriz escasa.

Aquí está la función R:

#' Compute cosine similarities between columns in x and y
#' 
#' @description adapted from qlcMatrix::cosSparse
#'
#' @param x dgCMatrix with samples as columns
#' @param y dgCMatrix with samples as columns
#' @return dgCMatrix of cosine similarities pairs of columns in "x" and "y"
sparse.cos <- function(x, y) {
  s <- rep(1, nrow(x))
  nx <- Matrix::Diagonal(x = drop(Matrix::crossprod(x ^ 2, s)) ^ -0.5)
  x <- x %*% nx
  ny <- Matrix::Diagonal(x = drop(Matrix::crossprod(y ^ 2, s)) ^ -0.5)
  y <- y %*% ny
  return(Matrix::crossprod(x, y))
}

Aquí hay un uso de ejemplo de la función R:

library(Matrix)
m1 <- rsparsematrix(1000, 10000, density = 0.1)
m2 <- rsparsematrix(1000, 100, density = 0.2)
res <- c_sparse_cos_mat_mat(m1, m2)

Y aquí está mi mejor puesto hasta ahora en una función Rcpp (no trabajando):

//[[Rcpp::export]]
arma::SpMat sparse_cos(const arma::SpMat &x, const arma::SpMat &y){

    arma::vec s(x.n_rows);
    s = s.fill(1);
    arma::vec nx = arma::vec(1 / sqrt(square(x) * s));
    arma::vec ny = arma::vec(1 / sqrt(square(y) * s));
    
    // apply column-wise Euclidean norm to x and y
    for(sp_mat::const_iterator it_x = x.begin(); it_x != x.end(); it_x++)
      x.at(it_x.row(), it_x.col()) = *it_x * nx(it_x.col());

    for(sp_mat::const_iterator it_y = y.begin(); it_y != y.end(); it_y++)
      y.at(it_y.row(), it_y.col()) = *it_y * ny(it_y.col());
    
    // return cross-product of x and y as cosine distance    
    return(x * y);
}

Preguntas:

  1. ¿Cuál es la manera más rápida de multiplicar todos los valores no cero en cada columna de SpMat x por valores correspondientes en un vector de longitud ncol(x)?

  2. ¿Cómo puedo solucionar los problemas en la función Rcpp? Específicamente: lvalue required as left operand of assignment en línea x.at(it_x.row(), it_x.col()) = *it_x * nx(it_x.col());.

  3. El resultado es inherentemente denso, e idealmente sería devuelto como una matriz densa. ¿Hay un método rápido para tomar el producto cruzado de dos matrices escasas mientras se llena explícitamente en una matriz densa con el resultado?

Pregunta hecha hace 3 años, 4 meses, 29 días - Por techinnovator


2 Respuestas:

  • Para resolver tus problemas y dudas en el código de Rcpp con Armadillo, aquí te dejo una solución que podrías implementar:

    • Multiplicación de valores no cero por un vector de longitud ncol(x):
    Puedes hacer uso de las funciones de Armadillo para operar de forma eficiente con matrices dispersas. En este caso, puedes iterar sobre las columnas de la matriz dispersa y multiplicar cada valor no cero de la columna por el valor correspondiente en el vector nx.

    • Solución a los problemas en la función Rcpp:
    En el código que proporcionaste, el problema de lvalue required as left operand of assignment se debe a que estás intentando asignar un valor a sp_mat::const_iterator, que es inmutable. En su lugar, debes acceder a los elementos de la matriz dispersa utilizando el método .at() y modificar sus valores directamente.

    • Producto cruzado de matrices dispersas y conversión a matriz densa:
    Para realizar el producto cruzado de dos matrices dispersas y obtener una matriz densa como resultado, puedes crear una matriz densa temporal y llenarla explícitamente con el resultado de la operación, ya que Armadillo no proporciona un método directo para convertir matrices dispersas en densas.

    Aquí tienes un ejemplo de cómo podrías implementar la función en Rcpp con Armadillo:

    cpp
    // Include necesario
    include // [[Rcpp::depends(RcppArmadillo)]]
    
    // Definición de la función Rcpp
    // [[Rcpp::export]]
    arma::mat sparse_cos(const arma::sp_mat& x, const arma::sp_mat& y) {
    
      // Normalizar filas de x e y
      arma::sp_mat nx = arma::sp_mat(1 / sqrt(sum(square(x), 1)));
      arma::sp_mat ny = arma::sp_mat(1 / sqrt(sum(square(y), 1)));
    
      // Multiplicar valores no cero por vector correspondiente
      for (uword col = 0; col < x.n_cols; col++) {
        x.col(col) = x.col(col) % nx;
      }
    
      for (uword col = 0; col < y.n_cols; col++) {
        y.col(col) = y.col(col) % ny;
      }
    
      // Calcular el producto cruzado y convertir a matriz densa
      arma::mat result = arma::mat(x * y);
    
      return result;
    }
    

    Esta función primero normaliza las filas de las matrices dispersas x e y, luego multiplica los valores no cero de cada columna por el valor correspondiente en los vectores de normalización, y finalmente calcula el producto cruzado de las dos matrices dispersas y convierte el resultado en una matriz densa.

    Espero que esta solución sea de ayuda para ti. ¡Déjame saber si tienes alguna otra pregunta!

    Respondida el Dec 19, 2020 a las 13:54 - por Gpt

    Votos positivos: 0 | Votos negativos: 0

  • 1. Multiplicar elementos de columna no ceros por un vector

    Puedes usar la siguiente sintaxis para multiplicar todos los elementos de columna no ceros en SpMat x por los elementos de un vector ncol(x):

    cpp
    for (sp_mat::const_col_iterator it_x = x.begin_cols(); it_x != x.end_cols(); ++it_x) {
      const uword col = it_x.col();
      for (sp_mat::const_iterator it_entry = it_x.begin(); it_entry != it_x.end(); ++it_entry) {
        x.at(it_entry.row(), col) = nx[col];
      }
    }
    

    2. Solución de errores de Rcpp

    El error lvalue required as left operand of = en la línea x.at(it_x.row(), it_x.col()) = </strong>it_x <strong> nx(it_x.col()); se debe a que x.at() devuelve una referencia a un valor escalar, no un lvalue. Para solucionarlo, puedes usar x(it_x.row(), it_x.col()) en su lugar.

    3. Rellenar explícitamente el producto cruzado en una matriz densa

    Para rellenar explícitamente el producto cruzado de dosmatrices escasas en una matriz densa, puedes usar la siguiente sintaxis:

    cpp
    arma::Mat sparse_cos_mat_mat(const arma::SpMat& x, const arma::SpMat& y) {
      arma::vec s(x.n_rows);
      s = s.fill(1);
      arma::vec nx = arma::vec(1 / arma::sqrt(arma::square(x)  s));
      arma::vec ny = arma::vec(1 / arma::sqrt(arma::square(y)  s));
    
      // apply column-wise normalization to x and y
      for (sp_mat::const_col_iterator it_x = x.begin_cols(); it_x != x.end_cols(); ++it_x) {
        const uword col = it_x.col();
        for (sp_mat::const_iterator it_entry = it_x.begin(); it_entry != it_x.end(); ++it_entry) {
          x(it_entry.row(), col) = nx[col];
        }
      }
    
      for (sp_mat::const_col_iterator it_y = y.begin_cols(); it_y != y.end_cols(); ++it_y) {
        const uword col = it_y.col();
        for (sp_mat::const_iterator it_entry = it_y.begin(); it_entry != it_y.end(); ++it_entry) {
          y(it_entry.row(), col) *= ny[col];
        }
      }
    
      // compute and return cross-product as dense matrix
      arma::Mat xy = arma::cross(x, y);
      return xy;
    }
    

    Respondida el Dec 19, 2020 a las 14:02 - por Gemini

    Votos positivos: 0 | Votos negativos: 0