Changing for loops to apply in r

Question

I have a nice random graphing simulation in a function that requires n nodes and a preferential attachment parameter beta. I use for loops, however when you take n to be very large, the code takes a long while to run. I was wondering if it was possible to use the apply family to make this more efficient.

binfunction <- function(y) { #Set up bins to create preferential attachment
 L <- length(y)
 x <- c(0, cumsum(y))
 U <- runif(1, min = 0 , max = sum(y))
  for(i in 1:L) {
   if(x[i] <= U && x[i+1] > U){
    return(i)
  }
 } 
}

random_graph <- function(n, beta) { #Random graphing function
 mat <- matrix(0,n,n)
 mat[1,2] <- 1
 mat[2,1] <- 1
  for(i in 3:n) {
   degvect <- colSums(mat[ , (1:(i-1))])
   degvect <- degvect^(beta)
   j <- binfunction(degvect)
   mat[i,j] <- 1
   mat[j,i] <- 1
 }
return(mat)
}

And it can be used with:

set.seed(123)
random_graph(10, 0.5)

could you please provide a [reproducible example](http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example/5963610#5963610) ? — Vincent Bonhomme, Apr 16 '16 at 04:55
@VincentBonhomme it is a random simulation. Just input `n` and `beta` into `random_graph`. (keep `beta` small, i.e. between 0 and 2) — Killian, Apr 16 '16 at 04:57
there is a 5% gain with just `j <- binfunction(colSums(mat[ , (1:(i-1))])^beta)`, yet not enough to solve your problem. — Vincent Bonhomme, Apr 16 '16 at 05:20
with the loop in `binfunction` replaced with `which(x > U)[1] - 1`, we have a 30% improvement (but not 100% sure it does the same job). — Vincent Bonhomme, Apr 16 '16 at 05:31
That is not exactly the same. Could you do `which(x > U & x < U)[1] - 1` ? Edit: that does not work either — Killian, Apr 16 '16 at 05:43
still ~35% faster. I use `microbenchmark(random_graph(100, 0.5), random_graph2(100, 0.5))` from `microbenchmark` package by the way. Don't think _I_ can do better but others may have clues. — Vincent Bonhomme, Apr 16 '16 at 05:47

score 0 · Answer 1 · answered Apr 16 '16 at 08:30

If you are looking for efficiency, vectorization is key.

The main time saver (at least for big n) is the use of sample.

For example:

n <- 1e7                             
sample(0:1, n, replace=TRUE)

takes about 0.2 seconds, while

for(i in 1:n) sample(0:1, 1)

takes about 24 sececonds. Vectorised operations can often replace loops, but knowing when and where depends on the familiarity with the available function for your needs.

score 0 · Answer 2 · answered Apr 16 '16 at 14:00

Consider using sparse matrix when dealing with large graph problem. Using normal matrix can be both computational and memory expensive. Here is a brief test by replacing your matrix with sparseMatrix from library Matrix:

 library(Matrix)
 psidom_random_graph <- function(n, beta) { #Random graphing function
     mat <- sparseMatrix(dims = c(n, n), i = {}, j = {})
     mat[1,2] <- T
     mat[2,1] <- T
     for(i in 3:n) {
         degvect <- colSums(mat[,1:(i-1)])
         degvect <- degvect^beta
         j <- binfunction(degvect)
         mat[i,j] <- T
         mat[j,i] <- T
    }
  return(mat)
}

> microbenchmark(psidom_random_graph(1000,1), killian_random_graph(1000,1))
Unit: seconds
                          expr      min       lq     mean   median       uq
  psidom_random_graph(1000, 1) 1.575727 1.593401 1.637430 1.602013 1.666144
 killian_random_graph(1000, 1) 7.031790 7.092761 7.212975 7.192756 7.308052
      max neval
 1.971281   100
 7.552074   100

Changing for loops to apply in r

2 Answers2