Skip to content

Is it possible to vectorize/accelerate this?

3 messages · hihi, michael.weylandt at gmail.com (R. Michael Weylandt, Dennis Murphy

hihi
# 
Dear Members,

I work on a simulaton experiment but it has an bottleneck. It's quite fast because of R and vectorizing, but it has a very slow for loop. The adjacent element of a vector (in terms of index number) depends conditionally on the former value of itself. Like a simple cumulating function (eg. cumsum) but with condition. Let's show me an example:
a_vec = rnorm(100)
b_vec = rep(0, 100)
b_vec[1]=a_vec[1]
for (i in 2:100){b_vec[i]=ifelse(abs(b_vec[i-1]+a_vec[i])>1, a_vec[i], b_vec[i-1]+a_vec[i])}
print(b_vec)

(The behaviour is like cumsum's, but when the value would excess 1.0 then it has another value from a_vec.)
Is it possible to make this faster? I experienced that my way is even slower than in Excel! Programming in C would my last try...
Any suggestions?

Than you,
Peter
michael.weylandt at gmail.com (R. Michael Weylandt
# 
I don't immediately see a good trick for vectorization so this seems to me to be a good candidate for work in a lower-level language. Staying within R, I'd suggest you use if and else rather than ifelse() since your computation isn't vectorized: this will eliminate a small amount over overhead. Since you also always add a_vec, you could also define b_vec as a copy of a to avoid all those calls to subset a, but I don't think the effects will be large and the code might not be as clear. 

You indicated that you may be comfortable with writing C, but I'd suggest you look into the Rcpp/Inline package pair which make the whole process much easier than it would otherwise be.

 I'm not at a computer write now or I'd write a fuller example, but the documentation for those packages is uncommonly good an you should be able to easily get it down into C++. If you aren't able to get it by tomorrow, let me know and I can help troubleshoot. The only things I foresee that you'll need to change are zero-basing, C's loop syntax, and (I think) the call to abs(). (I always forget where abs() lives in c++ ....)

The only possible hold up is that you need to be at a computer with a C compiler

Hope this helps,

Michael
On Nov 3, 2011, at 3:10 PM, hihi <v.p.mail at freemail.hu> wrote:

            

      
      
      
    

        

      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Dennis Murphy
    

    

      
      #
    

  


  

      
Hi:

You're doing the right thing in R by pre-allocating memory for the
result, but ifelse() is a vectorized function and your loop is
operating elementwise, so if-else is more appropriate. Try

for (i in 2:100){
    b_vec[i] <- if(abs(b_vec[i-1] + a_vec[i]) > 1) a_vec[i] else
b_vec[i-1] + a_vec[i]
  }

If speed is an issue, then I echo Michael's suggestion to write a
C(++) function and call it within R. The inline package is good for
this kind of thing.

HTH,
Dennis

        
On Thu, Nov 3, 2011 at 12:10 PM, hihi <v.p.mail at freemail.hu> wrote: