Ich würde es gerne als Zucker sehen. Leider bin ich nicht qualifiziert, es zu implementieren. Hier sind noch eine Reihe von verschiedenen Lösungen, mit denen ich gespielt habe.
Zuerst habe ich diese (colvec
statt vec
für tIdx
und Xmat.rows(...
statt X.rows(...
zur Arbeit kommen mussten einige Änderungen an Gong-Yi Liao Code machen:
mat Xmat(X.begin(), X.nrow(), X.ncol(), false);
colvec tIdx(T.begin(), T.size(), false);
mat y = Xmat.rows(find(tIdx == 1));
Zweitens sind hier drei Funktion mit Benchmarks, dass alle Subset-Matrizen auf einer logischen Anweisung basieren.Die Funktionen nehmen Argumente von arma oder rcpp an und geben Werte zurück Zwei basieren auf der Lösung von Gong-Yi Liao und eine ist eine einfache schleifenbasierte Lösung.
n (Zeilen) = 100, p (T == 1) = 0,3
expr min lq median uq max
1 submat_arma(X, T) 5.009 5.3955 5.8250 6.2250 28.320
2 submat_arma2(X, T) 4.859 5.2995 5.6895 6.1685 45.122
3 submat_rcpp(X, T) 5.831 6.3690 6.7465 7.3825 20.876
4 X[T == 1, ] 3.411 3.9380 4.1475 4.5345 27.981
n (Zeilen) = 10000, p (T == 1) = 0,3
expr min lq median uq max
1 submat_arma(X, T) 107.070 113.4000 125.5455 141.3700 1468.539
2 submat_arma2(X, T) 76.179 80.4295 88.2890 100.7525 1153.810
3 submat_rcpp(X, T) 244.242 247.3120 276.6385 309.2710 1934.126
4 X[T == 1, ] 229.884 236.1445 263.5240 289.2370 1876.980
submat.cpp
#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
using namespace Rcpp;
using namespace arma;
// arma in; arma out
// [[Rcpp::export]]
mat submat_arma(arma::mat X, arma::colvec T) {
mat y = X.rows(find(T == 1));
return y;
}
// rcpp in; arma out
// [[Rcpp::export]]
mat submat_arma2(NumericMatrix X, NumericVector T) {
mat Xmat(X.begin(), X.nrow(), X.ncol(), false);
colvec tIdx(T.begin(), T.size(), false);
mat y = Xmat.rows(find(tIdx == 1));
return y;
}
// rcpp in; rcpp out
// [[Rcpp::export]]
NumericMatrix submat_rcpp(NumericMatrix X, LogicalVector condition) {
int n=X.nrow(), k=X.ncol();
NumericMatrix out(sum(condition),k);
for (int i = 0, j = 0; i < n; i++) {
if(condition[i]) {
out(j,_) = X(i,_);
j = j+1;
}
}
return(out);
}
/*** R
library("microbenchmark")
# simulate data
n=100
p=0.3
T=rbinom(n,1,p)
X=as.matrix(cbind(rnorm(n),rnorm(n)))
# compare output
identical(X[T==1,],submat_arma(X,T))
identical(X[T==1,],submat_arma2(X,T))
identical(X[T==1,],submat_rcpp(X,T))
# benchmark
microbenchmark(X[T==1,],submat_arma(X,T),submat_arma2(X,T),submat_rcpp(X,T),times=500)
# increase n
n=10000
p=0.3
T=rbinom(n,1,p)
X=as.matrix(cbind(rnorm(n),rnorm(n)))
# benchmark
microbenchmark(X[T==1,],submat_arma(X,T),submat_arma2(X,T),submat_rcpp(X,T),times=500)
*/
ja, das Hinzufügen würde eine ganze Reihe von Entwicklung und Erprobung. Es ist also unwahrscheinlich, dass es bald passiert, es sei denn, es wird mit zweckgebundenen Mitteln finanziert –