returnOneR <- function(){
return(1)
}
#include <Rcpp.h>
int returnOneCpp(){
return(1);
}
int add(int, int, int);
int returnOneCpp();
> sourceCpp("returnOneCpp.cpp")
> returnOneCpp()
[1] 1
- Scalar input, scalar output
absR <- function(x){
if(x > 0){
x
}else if(x == 0){
0
}else{
-x
}
}
#include <Rcpp.h>
double absCpp(double x){
if(x > 0){
return x;
}else if (x == 0){
return 0;
}else{
return -x;
}
}
int add(int, int, int);
int returnOneCpp();
double absCpp(double x);
> sourceCpp("absCpp.cpp")
> absCpp(-3)
[1] 3
> absCpp(0)
[1] 0
> absCpp(4)
[1] 4
- Vector input, scalar output
- "Rcpp::NumericVector"という記述に注意
- ちなみに、cppFunction()を使うときには、double sumC(NumericVector x)と"Rcpp::"をつけなくてもうまくいきます…？？？
#include <Rcpp.h>
double meanCpp(Rcpp::NumericVector x){
int n = x.size();
double ret = 0;
for(int i = 0; i < n; ++i) {
ret += x[i];
}
return ret/n;
}
#include <Rcpp.h>
using namespace Rcpp;
double meanCpp2(NumericVector x){
int n = x.size();
double ret = 0;
for(int i = 0; i < n; ++i) {
ret += x[i];
}
return ret/n;
}
int add(int, int, int);
int returnOneCpp();
double absCpp(double);
double meanCpp(Rcpp::NumericVector);
meanR <- function(x){
n <- length(x)
ret <- 0
for(i in 1:length(x)){
ret <- ret + x[i]
}
ret/n
}
> sourceCpp("meanCpp.cpp")
> x <- runif(10)
> mean(x)
[1] 0.587014
> meanCpp(x)
[1] 0.587014
> meanR(x)
[1] 0.587014
> library(microbenchmark)
> x <- runif(1e3)
> microbenchmark(
+ mean(x),
+ meanR(x),
+ meanCpp(x)
+ )
Unit: microseconds
expr min lq median uq max neval
mean(x) 28.590 35.921 38.4870 48.567 110.328 100
meanR(x) 2275.825 2298.184 2389.4515 2869.430 3807.214 100
meanCpp(x) 7.331 10.080 12.4625 14.662 45.084 100
- Vector input, scalar output 2
- cppの関数の中からmeanCppを呼ぼうとしたけれど、うまくいかなかった。cppのmean,powの関数はmath.hをincludeしないといけないかと思ったが、不要らしい…
#include <Rcpp.h>
double corCpp(Rcpp::NumericVector x, Rcpp::NumericVector y) {
int n = x.size();
double mx = mean(x);
double my = mean(y);
double xx = 0;
double yy = 0;
double xy = 0;
for(int i = 0; i < n; ++i) {
xx += (x[i]-mx)*(x[i]-mx);
yy += (y[i]-my)*(y[i]-my);
xy += (x[i]-mx)*(y[i]-my);
}
return xy/pow(xx*yy,0.5);
}
int add(int, int, int);
int returnOneCpp();
double absCpp(double);
double meanCpp(Rcpp::NumericVector);
double corCpp(Rcpp::NumericVector,Rcpp::NumericVector);
corR <- function(x,y){
m.x <- mean(x)
m.y <- mean(y)
sum((x-m.x)*(y-m.y))/sqrt(sum((x-m.x)^2)*sum((y-m.y)^2))
}
> x <- runif(10)
> y <- runif(10)
> sourceCpp("corCpp.cpp")
> corCpp(x,y)
[1] -0.5522906
> corR(x,y)
[1] -0.5522906
> cor(x,y)
[1] -0.5522906
> x <- runif(1e3)
> y <- runif(1e3)
> microbenchmark(
+ cor(x,y),
+ corR(x,y),
+ corCpp(x,y)
+ )
Unit: microseconds
expr min lq median uq max neval
cor(x, y) 198.663 203.0620 212.2255 220.289 353.708 100
corR(x, y) 137.085 141.1170 147.8980 152.846 280.401 100
corCpp(x, y) 21.992 23.2755 27.8580 29.690 74.408 100
#include <Rcpp.h>
using namespace Rcpp;
NumericVector squareCpp(NumericVector x) {
int n = x.size();
NumericVector y(n);
for (int i = 0; i < n; i++){
y[i] = pow(x[i],2);
}
return y;
}
> sourceCpp("squareCpp2.cpp")
> x <- seq(from=0,to=1,length=4)
> x^2
[1] 0.0000000 0.1111111 0.4444444 1.0000000
> squareCpp(x)
[1] 0.0000000 0.1111111 0.4444444 1.0000000
#include <Rcpp.h>
using namespace Rcpp;
NumericVector sinCpp(NumericVector x) {
int n = x.size();
NumericVector y(n);
for (int i = 0; i < n; i++){
y[i] = sin(x[i]);
}
return y;
}
#include <Rcpp.h>
using namespace Rcpp;
NumericVector rowSumCpp(NumericMatrix x) {
int nr = x.nrow();
int nc = x.ncol();
NumericVector y(nr);
for (int i = 0; i < nr; i++){
double tmp = 0;
for(int j = 0; j < nc; j++){
tmp += x(i,j);
}
y[i] = tmp;
}
return y;
}
> x <- matrix(1:20,4,5)
> rowSumCpp(x)
[1] 45 50 55 60
> apply(x,1,sum)
[1] 45 50 55 60
> microbenchmark(
+ apply(x,1,sum),
+ rowSumCpp(x)
+ )
Unit: microseconds
expr min lq median uq max neval
apply(x, 1, sum) 1391.739 1433.891 1654.729 1775.320 5739.228 100
rowSumCpp(x) 40.686 42.336 50.216 52.049 103.364 100
>
