Jac*_*280 2 c++ r vector armadillo rcpp
我的问题涉及一个排序练习,我可以在R中轻松地(但可能很慢)进行,并希望用C++进行,以加快我的代码.
考虑三个相同大小的矢量a,b和c.在R中,以下命令首先按照b对数字进行排序,然后,在关系的情况下,将根据c进一步排序.
a<-a[order(b,c),1]
Run Code Online (Sandbox Code Playgroud)
例:
a<-c(1,2,3,4,5)
b<-c(1,2,1,2,1)
c<-c(5,4,3,2,1)
> a[order(b,c)]
[1] 5 3 1 4 2
Run Code Online (Sandbox Code Playgroud)
有没有一种有效的方法在C++中使用Armadillo向量进行此操作?
我们可以编写以下C++解决方案,我在一个文件中SO_answer.cpp:
#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
using namespace arma;
// [[Rcpp::export]]
vec arma_sort(vec x, vec y, vec z) {
// Order the elements of x by sorting y and z;
// we order by y unless there's a tie, then order by z.
// First create a vector of indices
uvec idx = regspace<uvec>(0, x.size() - 1);
// Then sort that vector by the values of y and z
std::sort(idx.begin(), idx.end(), [&](int i, int j){
if ( y[i] == y[j] ) {
return z[i] < z[j];
}
return y[i] < y[j];
});
// And return x in that order
return x(idx);
}
Run Code Online (Sandbox Code Playgroud)
我们所做的是利用std::sort()允许您根据自定义比较器进行排序的事实.我们使用比较器,z仅当元素y相等时才比较元素; 否则它比较的值y.1然后我们可以编译文件并在R中测试函数:
library(Rcpp)
sourceCpp("SO_answer.cpp")
set.seed(1234)
x <- sample(1:10)
y <- sample(1:10)
z <- sample(1:10)
y[sample(1:10, 1)] <- 1 # create a tie
all.equal(x[order(y, z)], c(arma_sort(x, y, z))) # check against R
# [1] TRUE # Good
Run Code Online (Sandbox Code Playgroud)
当然,我们还必须考虑这实际上是否会给你带来任何性能提升,这就是你做这件事的全部原因.我们的基准:
library(microbenchmark)
microbenchmark(r = x[order(y, z)],
arma = arma_sort(x, y, z),
times = 1e4)
Unit: microseconds
expr min lq mean median uq max neval cld
r 36.040 37.23 39.386160 37.64 38.32 3316.286 10000 b
arma 5.055 6.07 7.155676 7.00 7.53 107.230 10000 a
Run Code Online (Sandbox Code Playgroud)
在我的机器上,使用小向量看起来速度提高了5-6倍,但是当你向上扩展时这种优势并不好:
x <- sample(1:100)
y <- sample(1:100)
z <- sample(1:100)
y[sample(1:100, 10)] <- 1 # create some ties
all.equal(x[order(y, z)], c(arma_sort(x, y, z))) # check against R
# [1] TRUE # Good
microbenchmark(r = x[order(y, z)],
arma = arma_sort(x, y, z),
times = 1e4)
Unit: microseconds
expr min lq mean median uq max neval cld
r 44.50 46.360 48.01275 46.930 47.755 294.051 10000 b
arma 10.76 12.045 16.30033 13.015 13.715 5262.132 10000 a
x <- sample(1:1000)
y <- sample(1:1000)
z <- sample(1:1000)
y[sample(1:100, 10)] <- 1 # create some ties
all.equal(x[order(y, z)], c(arma_sort(x, y, z))) # check against R
# [1] TRUE # Good
microbenchmark(r = x[order(y, z)],
arma = arma_sort(x, y, z),
times = 1e4)
Unit: microseconds
expr min lq mean median uq max neval cld
r 113.765 118.7950 125.7387 120.5075 122.4475 3373.696 10000 b
arma 82.690 91.3925 104.0755 95.2350 99.4325 6040.162 10000 a
Run Code Online (Sandbox Code Playgroud)
它仍然更快,但是当你处于长度为1000的向量时,它的速度不到2倍.这可能就是为什么F.Privé说这个操作在R中应该足够快.虽然使用Rcpp转移到C++会给你带来很大的性能优势,你获得收益的程度在很大程度上取决于背景,正如Dirk Eddelbuettel在回答各种问题时多次提到的那样.
sort()或sort_index()(参见此处的Armadillo文档).如果您尝试vec按秒的值进行排序vec,可以x(arma::sort_index(y))按照我在此处的相关问题的答案中指出的那样使用.你甚至可以stable_sort_index()用来保护领带.但是,我无法弄清楚如何使用这些函数来解决您在此处提出的具体问题.