Rcpp Version 1.0.7
ifelseLooped.r
Go to the documentation of this file.
1 #!/usr/bin/env r
2 ##
3 ## This example goes back to the following StackOverflow questions:
4 ## http://stackoverflow.com/questions/7153586/can-i-vectorize-a-calculation-which-depends-on-previous-elements
5 ## and provides a nice example of how to accelerate path-dependent
6 ## loops which are harder to vectorise. It lead to the following blog
7 ## post:
8 ## http://dirk.eddelbuettel.com/blog/2011/08/23#rcpp_for_path_dependent_loops
9 ##
10 ## Thanks to Josh Ulrich for provided a first nice (R-based) answer on
11 ## StackOverflow and for also catching a small oversight in my posted answer.
12 ##
13 ## Dirk Eddelbuettel, 23 Aug 2011
14 ##
15 ## Copyrighted but of course GPL'ed
16
17
18 library(inline)
19 library(rbenchmark)
20 library(compiler)
21
22 fun1 <- function(z) {
23  for(i in 2:NROW(z)) {
24  z[i] <- ifelse(z[i-1]==1, 1, 0)
25  }
26  z
27 }
28 fun1c <- cmpfun(fun1)
29
30
31 fun2 <- function(z) {
32  for(i in 2:NROW(z)) {
33  z[i] <- if(z[i-1]==1) 1 else 0
34  }
35  z
36 }
37 fun2c <- cmpfun(fun2)
38
39
40 funRcpp <- cxxfunction(signature(zs="numeric"), plugin="Rcpp", body="
41  Rcpp::NumericVector z = Rcpp::NumericVector(zs);
42  int n = z.size();
43  for (int i=1; i<n; i++) {
44  z[i] = (z[i-1]==1.0 ? 1.0 : 0.0);
45  }
46  return(z);
47 ")
48
49
50 z <- rep(c(1,1,0,0,0,0), 100)
51 ## test all others against fun1 and make sure all are identical
52 all(sapply(list(fun2(z),fun1c(z),fun2c(z),funRcpp(z)), identical, fun1(z)))
53
54 res <- benchmark(fun1(z), fun2(z),
55  fun1c(z), fun2c(z),
56  funRcpp(z),
57  columns=c("test", "replications", "elapsed", "relative", "user.self", "sys.self"),
58  order="relative",
59  replications=1000)
60 print(res)
61
62 z <- c(1,1,0,0,0,0)
63 res2 <- benchmark(fun1(z), fun2(z),
64  fun1c(z), fun2c(z),
65  funRcpp(z),
66  columns=c("test", "replications", "elapsed", "relative", "user.self", "sys.self"),
67  order="relative",
68  replications=10000)
69 print(res2)
70
71
72 if (FALSE) { # quick test to see if Int vectors are faster: appears not
73  funRcppI <- cxxfunction(signature(zs="integer"), plugin="Rcpp", body="
74  Rcpp::IntegerVector z = Rcpp::IntegerVector(zs);
75  int n = z.size();
76  for (int i=1; i<n; i++) {
77  z[i] = (z[i-1]==1.0 ? 1.0 : 0.0);
78  }
79  return(z);
80  ")
81
82  z <- rep(c(1L,1L,0L,0L,0L,0L), 100)
83  identical(fun1(z),fun2(z),fun1c(z),fun2c(z),funRcppI(z))
84
85  res3 <- benchmark(fun1(z), fun2(z),
86  fun1c(z), fun2c(z),
87  funRcppI(z),
88  columns=c("test", "replications", "elapsed", "relative", "user.self", "sys.self"),
89  order="relative",
90  replications=1000)
91  print(res3)
92
93  z <- c(1L,1L,0L,0L,0L,0L)
94  res4 <- benchmark(fun1(z), fun2(z),
95  fun1c(z), fun2c(z),
96  funRcppI(z),
97  columns=c("test", "replications", "elapsed", "relative", "user.self", "sys.self"),
98  order="relative",
99  replications=10000)
100  print(res4)
101 }