Berechne euklidische Abstandsmatrix mit einem big.matrix Objekt

Question

Ich habe ein Objekt der Klasse big.matrix in R mit Dimension 778844 x 2. Die Werte sind alle Ganzzahlen (Kilometer). Mein Ziel ist es, die euklidische Distanzmatrix unter Verwendung der big.matrix zu berechnen und als Ergebnis ein Objekt der Klasse big.matrix zu haben. Ich würde gerne wissen, ob es einen optimalen Weg dafür gibt.

Der Grund für meine Wahl der Verwendung der Klasse big.matrix ist die Speicherbegrenzung. Ich könnte meine big.matrix in ein Objekt der Klasse matrix umwandeln und die euklidische Distanzmatrix mit dist() berechnen. dist() würde jedoch ein Objekt der Größe zurückgeben, die im Speicher nicht zugewiesen werden würde.

bearbeiten

Die folgende Antwort von John W. Emerson, Autor und Maintainer des bigmemory Paket gegeben wurde:

Sie große Algebra verwenden konnte ich erwarten, aber das wäre auch eine sein sehr schöner Anwendungsfall für Rcpp über sourceCpp(), und sehr kurz und einfach. Kurz gesagt, wir versuchen nicht einmal, High-Level-Features zu bieten (abgesehen von den Grundlagen, die wir als Proof-of-Concept implementiert haben). Kein einziger Algorithmus könnte alle Anwendungsfälle abdecken, wenn Sie beginnen, große Speicherprobleme zu vermeiden.

Answer 1

Hier ist ein Weg mit RcppArmadillo. Vieles davon ist dem RcppGallery example sehr ähnlich. Dies ergibt eine big.matrix mit den zugehörigen paarweisen (nach Zeilen) euklidischen Abständen. Ich mag meine big.matrix Funktionen in einer Wrapper-Funktion wickeln eine sauberere Syntax (dh vermeiden Sie die @address und andere Initialisierungen erstellen

. Hinweis - wie wir verwenden BigMemory (und daher besorgt mit RAM-Nutzung) Ich habe dieses Beispiel wieder die N-1 x N-1-Matrix von nur untere Dreieckselemente. Sie können dies ändern, aber das ist, was ich warf zusammen.

euc_dist.cpp

// To enable the functionality provided by Armadillo's various macros, 
// simply include them before you include the RcppArmadillo headers. 
#define ARMA_NO_DEBUG 

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

using namespace Rcpp; 
using namespace arma; 

// The following header file provides the definitions for the BigMatrix 
// object 
#include <bigmemory/BigMatrix.h> 

// C++11 plugin 
// [[Rcpp::plugins(cpp11)]] 

template <typename T> 
void BigArmaEuclidean(const Mat<T>& inBigMat, Mat<T> outBigMat) { 

    int W = inBigMat.n_rows; 

    for(int i = 0; i < W - 1; i++){ 
     for(int j=i+1; j < W; j++){ 
      outBigMat(j-1,i) = sqrt(sum(pow((inBigMat.row(i) - inBigMat.row(j)),2))); 
     } 
    } 
} 

// [[Rcpp::export]] 
void BigArmaEuc(SEXP pInBigMat, SEXP pOutBigMat) { 
    // First we tell Rcpp that the object we've been given is an external 
    // pointer. 
    XPtr<BigMatrix> xpMat(pInBigMat); 
    XPtr<BigMatrix> xpOutMat(pOutBigMat); 


    int type = xpMat->matrix_type(); 
    switch(type) { 
     case 1: 
     BigArmaEuclidean(
      arma::Mat<char>((char *)xpMat->matrix(), xpMat->nrow(), xpMat->ncol(), false), 
      arma::Mat<char>((char *)xpOutMat->matrix(), xpOutMat->nrow(), xpOutMat->ncol(), false) 
     ); 
     return; 

     case 2: 
     BigArmaEuclidean(
      arma::Mat<short>((short *)xpMat->matrix(), xpMat->nrow(), xpMat->ncol(), false), 
      arma::Mat<short>((short *)xpOutMat->matrix(), xpOutMat->nrow(), xpOutMat->ncol(), false) 
     ); 
     return; 

     case 4: 
     BigArmaEuclidean(
      arma::Mat<int>((int *)xpMat->matrix(), xpMat->nrow(), xpMat->ncol(), false), 
      arma::Mat<int>((int *)xpOutMat->matrix(), xpOutMat->nrow(), xpOutMat->ncol(), false) 
     ); 
     return; 

     case 8: 
     BigArmaEuclidean(
      arma::Mat<double>((double *)xpMat->matrix(), xpMat->nrow(), xpMat->ncol(), false), 
      arma::Mat<double>((double *)xpOutMat->matrix(), xpOutMat->nrow(), xpOutMat->ncol(), false) 
     ); 
     return; 

     default: 
     // We should never get here, but it resolves compiler warnings. 
     throw Rcpp::exception("Undefined type for provided big.matrix"); 
    } 

} 

Mein kleiner Wrapper

bigMatrixEuc <- function(bigMat){ 
    zeros <- big.matrix(nrow = nrow(bigMat)-1, 
         ncol = nrow(bigMat)-1, 
         init = 0, 
         type = typeof(bigMat)) 
    BigArmaEuc(bigMat@address, zeros@address) 
    return(zeros) 
} 

Der Test

library(Rcpp) 
sourceCpp("euc_dist.cpp") 

library(bigmemory) 

set.seed(123) 
mat <- matrix(rnorm(16), 4) 
bm <- as.big.matrix(mat) 

# Call new euclidean function 
bm_out <- bigMatrixEuc(bm)[] 

# pull out the matrix elements for out purposes 
distMat <- as.matrix(dist(mat)) 
distMat[upper.tri(distMat, diag=TRUE)] <- 0 
distMat <- distMat[2:4, 1:3] 

# check if identical 
all.equal(bm_out, distMat, check.attributes = FALSE) 
[1] TRUE