2015-02-10 8 views
10

में बार तर्क के साथ पुन: प्रजनन आर प्रतिनिधि मैं आरसीपीपी का उपयोग करना सीख रहा हूं। मैं आर 12+ में rep फ़ंक्शन को दोहराने के लिए सी ++ का उपयोग करना चाहता हूं। आर में कई rep के साथ कई चीनी फ़ंक्शन शामिल हैं (पृष्ठ 3 के नीचे देखें: http://cran.r-project.org/web/packages/Rcpp/vignettes/Rcpp-quickref.pdf। जैसा कि मैं प्रलेखन को समझता हूं, चीनी कार्य rep, rep_each, और rep_len दो तर्क लेते हैं - एक वेक्टर और एक पूर्णांक। हालांकि, मैं rep की कार्यक्षमता को दोहराता हूं जब मैं times तर्क का उपयोग करता हूं। उस स्थिति में, आप दो वैक्टरों की आपूर्ति कर सकते हैं। आर:सी ++ और आरसीपीपी

x <- c(10, 5, 12) 
y <- c(2, 6, 3) 

rep(x, times = y) 
[1] 10 10 5 5 5 5 5 5 12 12 12 

इस प्रकार times तर्क के साथ rep +०१२३९३१०१०८६ के प्रत्येक तत्व प्रतिकृतिसंबंधित y मान के रूप में कई बार। जैसा कि मैं इसे समझता हूं, मैं इसके लिए आरसीपीपी चीनी कार्यों का उपयोग करने का कोई तरीका नहीं देख सकता।

// [[Rcpp::export]] 
NumericVector reptest(NumericVector x, NumericVector y) { 
    int n = y.size(); 
    NumericVector myvector(sum(y)); 
    int ind = 0; 
    for (int i = 0; i < n; ++i) { 
     for (int j = 0; j < y(i); ++j) { 
      myvector(ind) = x[i]; 
      ind = ind + 1; 
      } 
     } 
    return myvector; 
} 

x <- c(10, 5, 12) 
y <- c(2, 6, 3) 

reptest(x, y) 
[1] 10 10 5 5 5 5 5 5 12 12 12 

यह थोड़ा आर में rep की तुलना में धीमी मैं अगर वहाँ वैसे भी है इस तेजी लाने के लिए या सोच रहा हूँ अगर किसी को भी एक बेहतर विचार है:

मैं निम्नलिखित सी ++ समारोह है कि काम करता है बनाया है। जैसा कि मैं इसे समझता हूं, rep सी कोड को कॉल कर रहा है, इसलिए शायद rep पर सुधार करना असंभव होगा। मेरा लक्ष्य एक एमसीएमसी पाश को गति देना है (जो rep फ़ंक्शन का उपयोग करता है) जो आर में चलाने के लिए बहुत समय लगता है, इसलिए कोई भी गति उपयोगी होगी। एमसीएमसी पाश के अन्य भाग धीमे भाग हैं, rep नहीं, लेकिन मुझे अपने लूप में एक ही कार्यक्षमता की आवश्यकता है। यह इसकी गति बढ़ाने के लिए

+1

सुनिश्चित नहीं हैं कि अगर उपयोगी है, लेकिन यहां 'rep' स्रोत कोड के लिए एक कड़ी है। ऐसा लगता है कि यह 'सी' में है। https://github.com/wch/r-source/blob/ed415a8431b32e079100f50a846e4769aeb54d5a/src/main/seq.c –

+1

बहुत * बहुत * पहली त्वरित नज़र, यह ठीक दिखता है। आप अधिक यथार्थवादी आकारों पर बेंचमार्क करना चाहते हैं। –

उत्तर

9

यहां दो आंतरिक संस्करणों पर एक त्वरित रैफ है। यह भी rep.int() जोड़ता है:

#include <algorithm> 
#include <Rcpp.h> 
using namespace Rcpp; 

// [[Rcpp::export]] 
NumericVector reptest(NumericVector x, NumericVector y) { 
    int n = y.size(); 
    NumericVector myvector(sum(y)); 
    int ind = 0; 
    for (int i = 0; i < n; ++i) { 
     for (int j = 0; j < y[i]; ++j) { 
      myvector[ind] = x[i]; 
      ind = ind + 1; 
      } 
     } 
    return myvector; 
} 

// [[Rcpp::export]] 
NumericVector reptest2(NumericVector x, NumericVector y) { 
    int n = y.size(); 
    NumericVector myvector(sum(y)); 
    int ind=0; 
    for (int i=0; i < n; ++i) { 
     int p = y[i]; 
     std::fill(myvector.begin()+ind, myvector.begin()+ind+p, x[i]); 
     ind += p; 
    } 
    return myvector; 
} 


/*** R 
x <- rep(c(10, 5, 12), 10000) 
y <- rep(c(20, 60, 30), 10000) 
all.equal(reptest(x, y), reptest2(x, y), rep(x, times=y)) 

library(microbenchmark) 
microbenchmark(reptest(x, y), reptest2(x, y), rep(x, times=y), rep.int(x, y)) 
***/ 

इस के साथ, हम एक छोटे से करीब लेकिन आर अभी भी जीत:

// [[Rcpp::plugins(cpp11)]] 
// [[Rcpp::export]] 
NumericVector reptest3(const NumericVector& x, const IntegerVector& times) { 
    std::size_t n = times.size(); 
    if (n != 1 && n != x.size()) 
     stop("Invalid 'times' value"); 
    std::size_t n_out = std::accumulate(times.begin(), times.end(), 0); 
    NumericVector res = no_init(n_out); 
    auto begin = res.begin(); 
    for (std::size_t i = 0, ind = 0; i < n; ind += times[i], ++i) { 
     auto start = begin + ind; 
     auto end = start + times[i]; 
     std::fill(start, end, x[i]); 
    } 
    return res; 
} 

:

R> Rcpp::sourceCpp("/tmp/rep.cpp") 

R> x <- rep(c(10, 5, 12), 10000) 

R> y <- rep(c(20, 60, 30), 10000) 

R> all.equal(reptest(x, y), reptest2(x, y), rep(x, times=y)) 
[1] TRUE 

R> library(microbenchmark) 

R> microbenchmark(reptest(x, y), reptest2(x, y), rep(x, times=y), rep.int(x, y)) 
Unit: milliseconds 
       expr  min  lq mean median  uq  max neval 
    reptest(x, y) 4.61604 4.74203 5.47543 4.78120 6.78039 7.01879 100 
    reptest2(x, y) 3.14788 3.27507 5.25515 3.33166 5.24583 140.64080 100 
rep(x, times = y) 2.45876 2.56025 3.26857 2.60669 4.60116 6.76278 100 
    rep.int(x, y) 2.42390 2.50241 3.38362 2.53987 4.56338 6.44241 100 
R> 
9

एक तरीका यह बजाय प्रत्येक तत्व के माध्यम से पुनरावृत्ति भरे जाने की std::fill उपयोग करने के लिए होगा:

#include <algorithm> 
#include <Rcpp.h> 
using namespace Rcpp; 
// [[Rcpp::export]] 
NumericVector reptest2(NumericVector x, NumericVector y) { 
    int n = y.size(); 
    std::vector<double> myvector(sum(y)); 
    int ind=0; 
    for (int i=0; i < n; ++i) { 
    std::fill(myvector.begin()+ind, myvector.begin()+ind+y[i], x[i]); 
    ind += y[i]; 
    } 
    return Rcpp::wrap(myvector); 
} 

एक बड़ा उदाहरण पर, यह rep के साथ नज़दीकी बढ़ाने प्रकट होता है:

x <- rep(c(10, 5, 12), 10000) 
y <- rep(c(20, 60, 30), 10000) 
all.equal(reptest(x, y), reptest2(x, y), rep(x, times=y)) 
# [1] TRUE 

library(microbenchmark) 
microbenchmark(reptest(x, y), reptest2(x, y), rep(x, times=y)) 
# Unit: milliseconds 
#    expr  min  lq  mean median  uq  max neval 
#  reptest(x, y) 9.072083 9.297573 11.469345 9.522182 13.015692 20.47905 100 
#  reptest2(x, y) 5.097358 5.270827 7.367577 5.436549 8.961004 15.68812 100 
# rep(x, times = y) 1.457933 1.499051 2.884887 1.561408 1.949750 13.21706 100 
+0

अच्छा, इस बारे में भी सोचा। हालांकि आपको अंतिम 'रैप() 'की आवश्यकता नहीं है। या इसके बजाय, यदि आप एक प्रतिलिपि भी मजबूर करते हैं। –

+0

'my_vector'' std :: vector' क्यों है? इसकी याददाश्त को 'न्यूमेरिक वेक्टर' में कॉपी करना होगा। पहली जगह में 'न्यूमेरिक वेक्टर' से शुरू क्यों न करें। –

4

हम no_init साथ आर आधार rep प्रदर्शन प्राप्त कर सकते हैं बेंचमार्क:

library(microbenchmark) 
x <- rep(c(10, 5, 12), 10000) 
y <- rep(c(20, 60, 30), 10000) 
microbenchmark(
    reptest(x, y), reptest2(x, y), reptest3(x, y), 
    rep(x, times = y), rep.int(x, y)) 
#> Unit: milliseconds 
#>    expr  min  lq  mean median  uq  max neval 
#>  reptest(x, y) 13.209912 14.014886 15.129395 14.457418 15.123676 56.655527 100 
#> reptest2(x, y) 4.289786 4.653088 5.789094 5.105859 5.782284 46.679824 100 
#> reptest3(x, y) 1.812713 2.810637 3.860590 3.194529 3.809141 44.111422 100 
#> rep(x, times = y) 2.510219 2.877324 3.576183 3.461315 3.927312 5.961317 100 
#>  rep.int(x, y) 2.496481 2.901303 3.422384 3.318761 3.831794 5.283187 100 

इसके अलावा हम RcppParallel के साथ इस कोड में सुधार कर सकते हैं:

struct Sum : Worker { 
    const RVector<int> input; 
    int value; 
    Sum(const IntegerVector& input) : input(input), value(0) {} 
    Sum(const Sum& sum, Split) : input(sum.input), value(0) {} 
    void operator()(std::size_t begin, std::size_t end) { 
     value += std::accumulate(input.begin() + begin, input.begin() + end, 0); 
    } 
    void join(const Sum& rhs) { 
     value += rhs.value; 
    } 
}; 

struct Fill: Worker { 
    const RVector<double> input; 
    const RVector<int> times; 
    RVector<double> output; 
    std::size_t ind; 
    Fill(const NumericVector& input, const IntegerVector& times, NumericVector& output) 
     : input(input), times(times), output(output), ind(0) {} 
    void operator()(std::size_t begin, std::size_t end) { 
     for (std::size_t i = begin; i < end; ind += times[i], ++i) 
      std::fill(output.begin() + ind, output.begin() + ind + times[i], input[i]); 
    } 
}; 

// [[Rcpp::export]] 
NumericVector reptest4(const NumericVector& x, const IntegerVector& times) { 
    std::size_t n = times.size(); 
    if (n != 1 && n != x.size()) 
     stop("Invalid 'times' value"); 
    Sum s(times); 
    parallelReduce(0, n, s); 
    NumericVector res = no_init(s.value); 
    Fill f(x, times, res); 
    parallelFor(0, n, f); 
    return res; 
} 

तुलना:

library(microbenchmark) 
x <- rep(c(10, 5, 12), 10000) 
y <- rep(c(20, 60, 30), 10000) 
microbenchmark(
    reptest(x, y), reptest2(x, y), reptest3(x, y), 
    rep(x, times = y), rep.int(x, y)) 
#> Unit: milliseconds 
#>    expr  min  lq  mean median  uq  max neval 
#>  reptest3(x, y) 2.442446 3.410985 5.143627 3.893345 5.054285 57.871429 100 
#>  reptest4(x, y) 1.211256 1.534428 1.979526 1.821398 2.170999 4.073395 100 
#> rep(x, times = y) 2.435122 3.173904 4.447954 3.795285 4.687695 54.000920 100 
#>  rep.int(x, y) 2.444310 3.208522 4.026722 3.913618 4.798793 6.690333 100 
संबंधित मुद्दे