R schnelle Cbind-Matrix mit Rcpp

Question

cbind in R ist relativ zeitaufwendig in wiederholten Aufrufen, aber es ist auch für verschiedene Datentypen mächtig. Ich habe Code geschrieben, der 3x schneller ist als cbind beim Binden von zwei Matrizen. Aber bind_cols in dplyr Paket ist nur 100X schneller als cbind. Es ist nur schade, dass es Matrix nicht als Eingabe nehmen kann. Kann jemand den folgenden Code schneller machen? Auch, wie binde ich schnell spärliche Matrix? Hier ist der Code, den ich verwendet:

require(Rcpp) 

func <- 'NumericMatrix mmult(NumericMatrix a,NumericMatrix b) { 
    //the colnumber of first matrix 
    int acoln=a.ncol(); 
    //the colnumber of second matrix 
    int bcoln=b.ncol(); 
    //build a new matrix, the dim is a.nrow() and acoln+bcoln 
    NumericMatrix out(a.nrow(),acoln+bcoln) ; 
    for (int j = 0; j < acoln + bcoln; j++) { 
     if (j < acoln) { 
      out(_,j) = a(_,j); 
     } else { 
      //put the context in the second matrix to the new matrix 
      out(_,j) = b(_,j-acoln); 
     } 
    } 
    return out ; 
}' 

a <- matrix(rep(1,2000*100),2000) 
b <- matrix(rep(2,2000*10),2000) 

cppFunction(func) 

system.time(for (i in seq(1,800)) {mmult(a,b)}) 
system.time(for (i in seq(1,800)) {cbind(a,b)}) 
identical(mmult(a,b),cbind(a,b)) 

Answer 1

Borrowing eine Idee von diesem comment von Romain Francois auf einem meiner früheren RCPP Abenteuer,

func1 <- 'NumericMatrix mmult1(NumericMatrix a, NumericMatrix b) { 
    int acoln = a.ncol(); 
    int bcoln = b.ncol(); 
    NumericMatrix out = no_init_matrix(a.nrow(), acoln + bcoln); 
    for (int j = 0; j < acoln + bcoln; j++) { 
    if (j < acoln) { 
     out(_, j) = a(_, j); 
    } else { 
     out(_, j) = b(_, j - acoln); 
    } 
    } 
    return out; 
}' 

cppFunction(func1) 
set.seed(42) 
a <- matrix(rnorm(1e7), 1e3) 
b <- matrix(runif(1e7), 1e3) 

identical(mmult(a, b), mmult1(a, b)) 
#TRUE 

library(microbenchmark) 
microbenchmark(mmult(a, b), 
       mmult1(a, b), 
       cbind(a, b), 
       times = 10) 

#Unit: milliseconds 
#   expr min  lq mean median uq max neval 
# mmult(a, b) 69.64 70.52 89.71 72.28 128.8 136.6 10 
# mmult1(a, b) 50.84 50.95 69.65 51.43 111.6 114.4 10 
# cbind(a, b) 192.35 194.67 201.13 195.30 196.1 255.9 10 

Keine große Sache, aber auch nicht schlecht für eine solche triviale Veränderung.