Подмножество Rcpp Matrix, которое соответствует логическому утверждению

В R, если у нас есть матрица данных, скажем, матрица X 100 × 10, и вектор t из 100 элементов с возможными значениями (0, 1, 2, 3), мы можем легко найти подматрицу y из X, используя простую синтаксис:

y = X[t == 1, ]

Но проблема в том, как я могу сделать это с NumericMatrix Rcpp?
(Или, в более общем смысле, как я могу сделать это в любых контейнерах C ++?)

Благодаря подсказке Дирка кажется, что

NumericMatrix X(dataX);
IntegerVector T(dataT);
mat Xmat(X.begin(), X.nrow(), X.ncol(), false);
vec tIdx(T.begin(), T.size(), false);
mat y = X.rows(find(tIdx == 1));

Можно делать то, что я хочу, но это кажется слишком длинным.

16

Решение

Наиболее близким из известных мне является сочетание find() функция в сочетании с submat() функция в броненосец доступны через RcppArmadillo.

Редактировать: Это, конечно, то, что мы могли бы добавить через патч. Если кто-то достаточно мотивирован, чтобы попробовать это, пожалуйста, зайдите в список рассылки rcpp-devel.

8

Другие решения

Я хотел бы видеть это как сахар. К сожалению, я не квалифицирован, чтобы реализовать это все же. Вот еще ряд различных решений, с которыми я играл.

Сначала мне пришлось внести некоторые изменения в код Gong-Yi Liao, чтобы заставить это работать (colvec вместо vec за tIdx а также Xmat.rows(... вместо X.rows(...:

mat Xmat(X.begin(), X.nrow(), X.ncol(), false);
colvec tIdx(T.begin(), T.size(), false);
mat y = Xmat.rows(find(tIdx == 1));

Во-вторых, здесь есть три функции с эталонами, которые все подмножества матриц основаны на логическом утверждении. Функции принимают аргументы arma или rcpp и возвращают значения Два основаны на решении Гонг-И Ляо, а одна — простое решение на основе цикла.

n (ряды) = 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 (ряды) = 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)

*/
8

По вопросам рекламы [email protected]