R help - Nov 2011 - Is it possible to vectorize/accelerate this?

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

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:
> 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
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

An embedded and charset-unspecified text was scrubbed...
Name: nem el?rhet?
URL:
<https://stat.ethz.ch/pipermail/r-help/attachments/20111103/536f0c1a/attachment.pl>

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:
> 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
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

Yes -- if & else is much faster than ifelse() because if is a
primitive while ifelse() is a whole function call (in fact, you can
see the code by typing ifelse into the prompt and see that it has two
if calls within it.

Michael

On Thu, Nov 3, 2011 at 4:38 PM, hihi <v.p.mail at freemail.hu>
wrote:
> Hi,
> thank you for your very immediate response. :-) Is if than and else faster
> than ifelse? I'm wondering (or not knowing something)
> Best regards,
> Peter
> 2011/11/3 R. Michael Weylandt <michael.weylandt at gmail.com>
> <michael.weylandt at gmail.com>
>>
>> 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:
>>
>> > 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
>> >
>> > ______________________________________________
>> > R-help at r-project.org mailing list
>> > https://stat.ethz.ch/mailman/listinfo/r-help
>> > PLEASE do read the posting guide
>> > http://www.R-project.org/posting-guide.html
>> > and provide commented, minimal, self-contained, reproducible code.
>
>

I have to admit to not doing careful timing tests, but I often eliminate 
if() lines as follows (bad/good is just my preference)

BAD:  b[i] <- if(a[i]>1) a[i] else a[i-1]

GOOD:  b[i] <- a[i]* (a[i]>1) + a[i-1] * (a[i]<=1)




On Thu, Nov 3, 2011 at 12:10 PM, hihi <v.p.mail_at_freemail.hu> wrote:
 > 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)
-- 

Sent from my Cray XK6
"Pendeo-navem mei anguillae plena est."