[R] Is R the right choice for simulating first passage times of random walks?

[Steve Lianoglou]

Hi,

I haven't been following this thread very closely, but I'm getting the
impression that the "inner loop" that's killing you folks here looks
quite simple (assuming it is the one provided below).

How about trying to write the of this `f4` function below using the
rcpp/inline combo. The C/C++ you will need to write looks to be quite
trivial, let's change f4 to accept an x argument as a vector:

I've defined f4 in the same way as Dennis did:

> f4 <- function()
>  {
>     x <- sample(c(-1L,1L), 1)
>
>      if (x >= 0 ) {return(1)} else {
>           csum <- x
>           len <- 1
>               while(csum < 0) {
>                   csum <- csum + sample(c(-1, 1), 1)
>                   len <- len + 1
>                                          }     }
>      len
>  }

Now, let's do some inline/c++ mojo:

library(inline)
inc <- "
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
"

fxx <-cxxfunction(includes=inc, plugin="Rcpp", body="
  int len = 1;
  int x = ((rand() % 2 ) == 0) ? 1 : -1;
  int csum = x;

  while (csum < 0) {
    x = ((rand() % 2 ) == 0) ? 1 : -1;
    len++;
    csum = csum + x;
  }

  return wrap(len);
")

Assuming I've faithfully translated this into c++, the timings aren't
all that comparable.

Doing 500 replicates with the pure R version:

set.seed(123)
system.time(out <- replicate(500, f4()))
   user  system elapsed
 31.525   0.120  32.510

Doing 10,000 replicates using the fxx function doesn't even break a sweat:

system.time(outxx <- replicate(10000, fxx()))
   user  system elapsed
  0.371   0.001   0.373

range(out)
[1]       1 1994308

range(outxx)
[1]        1 11909394

Hope I'm not too off of the mark, here.
-steve
-- 
Steve Lianoglou
Graduate Student: Computational Systems Biology
 | Memorial Sloan-Kettering Cancer Center
 | Weill Medical College of Cornell University
Contact Info: http://cbio.mskcc.org/~lianos/contact