2014-01-26 15 views
15

是否有任何C++函數等效於R中%in%操作?考慮下面的命令在R:A C++的%版本%運算中的R

which(y %in% x) 

我試圖找到用C東西相當於++(特別是在犰狳),我找不到任何東西。然後我寫了自己的函數,與上面的R命令相比,它非常慢。

這裏是我寫的:

#include <RcppArmadillo.h> 
// [[Rcpp::depends("RcppArmadillo")]] 

// [[Rcpp::export]] 
arma::uvec myInOperator(arma::vec myBigVec, arma::vec mySmallVec){ 
arma::uvec rslt = find(myBigVec == mySmallVec[0]); 
for (int i = 1; i < mySmallVec.size(); i++){ 
    arma::uvec rslt_tmp = find(myBigVec == mySmallVec[i]); 
    rslt = arma::unique(join_cols(rslt, rslt_tmp)); 
} 
return rslt; 
} 

現在在上面的代碼採購後,我們有:

x <- 1:4 
y <- 1:10 
res <- benchmark(myInOperator(y, x), which(y %in% x), columns = c("test", 
     "replications", "elapsed", "relative", "user.self", "sys.self"), 
     order = "relative") 

而且這裏的結果:

    test replications elapsed relative user.self sys.self 
2 which(y %in% x)   100 0.001  1  0.001  0 
1 myInOperator(y, x)   100 0.002  2  0.001  0 

任何人都可以指導我尋找一個相應的C++代碼(y%in%x)或讓我的代碼更高效?這兩種功能的使用時間已經非常短了。我想我的效率意味着更多的是從編程的角度,以及我對問題的看法以及我使用的命令是否有效。

我感謝您的幫助。

+1

這將是硬的(在%x和y%)打''哪,因爲' %in(調用'match')和'which'已經是'.Internal'函數,因此在C中實現並可能被優化。雖然可以通過避免'%in%'生成的臨時邏輯向量來提高性能。 –

+1

它可能會有點幫助,如果你通過引用而不是複印,'(常量ARMA :: VEC&myBigVec,常量ARMA :: VEC&mySmallVec)傳遞的東西' – OneOfOne

+2

我dont't不夠了解RCPP(並沒有什麼的阿米迪洛),所以我不能回答這個問題。但是,如果我在C++中這樣做,我會查找'std :: set_intersection'。 –

回答

15

編輯:感謝@MatthewLundberg和@Yakk追趕我的愚蠢錯誤。

如果你真正想要的只是更快的匹配,你應該檢查出西門Urbanek的fastmatch包。然而,Rcpp確實有一個糖in功能,可以在這裏使用。 in使用了fastmatch包中的一些想法,並將它們合併到Rcpp。我也在這裏比較@ hadley的解決方案。

// [[Rcpp::plugins("cpp11")]] 
#include <Rcpp.h> 
using namespace Rcpp; 

// [[Rcpp::export]] 
std::vector<int> sugar_in(IntegerVector x, IntegerVector y) { 
    LogicalVector ind = in(x, y); 
    int n = ind.size(); 
    std::vector<int> output; 
    output.reserve(n); 
    for (int i=0; i < n; ++i) { 
    if (ind[i]) output.push_back(i+1); 
    } 
    return output; 
} 

// [[Rcpp::export]] 
std::vector<int> which_in(IntegerVector x, IntegerVector y) { 
    int nx = x.size(); 
    std::unordered_set<int> z(y.begin(), y.end()); 
    std::vector<int> output; 
    output.reserve(nx); 
    for (int i=0; i < nx; ++i) { 
    if (z.find(x[i]) != z.end()) { 
     output.push_back(i+1); 
    } 
    } 
    return output; 
} 


// [[Rcpp::export]] 
std::vector<int> which_in2(IntegerVector x, IntegerVector y) { 
    std::vector<int> y_sort(y.size()); 
    std::partial_sort_copy (y.begin(), y.end(), y_sort.begin(), y_sort.end()); 

    int nx = x.size(); 
    std::vector<int> out; 

    for (int i = 0; i < nx; ++i) { 
    std::vector<int>::iterator found = 
     lower_bound(y_sort.begin(), y_sort.end(), x[i]); 
    if (found != y_sort.end()) { 
     out.push_back(i + 1); 
    } 
    } 
    return out; 
} 

/*** R 
set.seed(123) 
library(microbenchmark) 
x <- sample(1:100) 
y <- sample(1:10000, 1000) 
identical(sugar_in(y, x), which(y %in% x)) 
identical(which_in(y, x), which(y %in% x)) 
identical(which_in2(y, x), which(y %in% x)) 
microbenchmark(
    sugar_in(y, x), 
    which_in(y, x), 
    which_in2(y, x), 
    which(y %in% x) 
) 
*/ 

調用此sourceCpp給我,從基準,

Unit: microseconds 
      expr min  lq median  uq max neval 
    sugar_in(y, x) 7.590 10.0795 11.4825 14.3630 32.753 100 
    which_in(y, x) 40.757 42.4460 43.4400 46.8240 63.690 100 
which_in2(y, x) 14.325 15.2365 16.7005 17.2620 30.580 100 
which(y %in% x) 17.070 21.6145 23.7070 29.0105 78.009 100 
+1

當你調用'which_in'時,你必須複製矢量嗎? –

+1

lhs使得返回值順序錯誤。嘗試使用從rhs複製的無序集,並轉換lhs而不先複製它?更快的對數因子。額外的信貸產生懶惰的輸出(mayhap使用'boost')。哦,並注意任何可迭代的範圍是一個有效的lhs和rhs:不知道如果R使它無用。 – Yakk

+1

根據[此Rcpp圖庫文章](http://gallery.rcpp.org/articles/armadillo-subsetting/),Armadillo的find()如何? –

9

對於這組輸入,我們可以通過使用在技術上具有較高的算法複雜度的方法勉強維持多一點表現( O(ln n)與O(1)進行查找),但具有較低的常量:二進制搜索。

// [[Rcpp::plugins("cpp11")]] 
#include <Rcpp.h> 
using namespace Rcpp; 

// [[Rcpp::export]] 
std::vector<int> which_in(IntegerVector x, IntegerVector y) { 
    int nx = x.size(); 
    std::unordered_set<int> z(y.begin(), y.end()); 
    std::vector<int> output; 
    output.reserve(nx); 
    for (int i=0; i < nx; ++i) { 
    if (z.find(x[i]) != z.end()) { 
     output.push_back(i+1); 
    } 
    } 
    return output; 
} 

// [[Rcpp::export]] 
std::vector<int> which_in2(IntegerVector x, IntegerVector y) { 
    std::vector<int> y_sort(y.size()); 
    std::partial_sort_copy (y.begin(), y.end(), y_sort.begin(), y_sort.end()); 

    int nx = x.size(); 
    std::vector<int> out; 

    for (int i = 0; i < nx; ++i) { 
    std::vector<int>::iterator found = 
     lower_bound(y_sort.begin(), y_sort.end(), x[i]); 
    if (found != y_sort.end()) { 
     out.push_back(i + 1); 
    } 
    } 
    return out; 
} 

/*** R 
set.seed(123) 
library(microbenchmark) 
x <- sample(1:100) 
y <- sample(1:10000, 1000) 
identical(which_in(y, x), which(y %in% x)) 
identical(which_in2(y, x), which(y %in% x)) 
microbenchmark(
    which_in(y, x), 
    which_in2(y, x), 
    which(y %in% x) 
) 
*/ 

在我的電腦能產生

Unit: microseconds 
      expr min lq median uq max neval 
    which_in(y, x) 39.3 41.0 42.7 44.0 81.5 100 
which_in2(y, x) 12.8 13.6 14.4 15.0 23.8 100 
which(y %in% x) 16.8 20.2 21.0 21.9 31.1 100 

所以比基地R.更好的約30%

+1

你的意思是說'which_in2'具有複雜性O(n log n)? –

+0

@MatthewLundberg哎呀,我感到困惑。我的意思是O(1),指的是個人查詢。 – hadley