Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Avoiding for loop in R

I'd like some advice on writing better R code. I have currently written a loop in R and it suffers from performance issues.

I can't wrap my head around vectorizing it, because each row in the output data frame has dependencies upon earlier rows and they trickle down iteratively, so I have written a loop to read/write the rows in sequence.

An example of my code:

example <- data.frame(a=c(.5,.1,.5,.25),b=c(1,0,2,0),c=c(1,2,3,4),d=c(4,3,2,1))

for (i in 2:nrow(example)) {
  if (example[i,1]>0) {
    example[i,2]<-example[i,2]+example[i-1,2]*example[i,1]
    example[i,3]<-example[i,3]+example[i-1,3]*example[i,1]
    example[i,4]<-example[i,4]+example[i-1,4]*example[i,1]

  }
}

To see what's happening:

# before    
     a b c d
1 0.50 1 1 4
2 0.10 0 2 3
3 0.50 2 3 2
4 0.25 0 4 1

# after
     a      b      c     d
1 0.50 1.0000 1.0000 4.000
2 0.10 0.1000 2.1000 3.400
3 0.50 2.0500 4.0500 3.700
4 0.25 0.5125 5.0125 1.925
like image 268
adamrussell Avatar asked Aug 02 '16 16:08

adamrussell


1 Answers

I'm not sure how to avoid by row operations, but here are 3 advises that will improve performance by ~X90

  1. Use matrices instead of data.frames for such operations
  2. Pre-alocate the rows you want to operate on
  3. Vectorize operations within the loop

In other words, try converting your code to

indx <- which(example[-1, 1] > 0)
for(i in indx + 1) example[i, -1] <- example[i, -1] + example[i-1, -1] * example[i, 1]

Also, note that this solution is generalized for any number of columns


Benchmark

set.seed(123)
N <- 1e3
test <- matrix(runif(N * 4), ncol = 4)
example <- as.data.frame(test)


OP <- function(x){
  for (i in 2:nrow(x)) {
    if (x[i, 1]>0) {
      x[i,2]<-x[i,2]+x[i-1,2]*x[i,1]
      x[i,3]<-x[i,3]+x[i-1,3]*x[i,1]
      x[i,4]<-x[i,4]+x[i-1,4]*x[i,1]
    }
  }
  x
}

David <- function(x){
  indx <- which(x[-1, 1] > 0)
  for(i in indx + 1) x[i, -1] <- x[i, -1] + x[i-1, -1] * x[i, 1]
  x
}

identical(OP(example), as.data.frame(David(test)))
# [1] TRUE

library(microbenchmark)    
microbenchmark(OP(example), David(test))    
# Unit: milliseconds
#        expr        min         lq       mean     median         uq        max neval cld
# OP(example) 243.913429 246.248061 257.672703 247.104350 256.701590 337.375850   100   b
# David(test)   3.020688   3.080685   3.336778   3.133483   3.301797   9.240615   100  a 
like image 60
David Arenburg Avatar answered Nov 13 '22 11:11

David Arenburg