How to optimize my R script in order to use “multicore”

Question

I am using GNU R at a Ubuntu-Lucid PC which has 4 CPUs. In order to use all 4 CPUs, I installed the "r-cran-multicore" package. As the package's manual lacks of practical examples that I understand, I need advice in how to optimize my script in order to make use of all 4 CPUs.

My dataset is a data.frame (called P1) that has 50,000 rows and 1600 cols. For each row, I'd like to calc the maximun, sum and mean. My script looks as follows:

p1max <- 0
p1mean <- 0
p1sum <-0
plength <- length(P1[,1])
for(i in 1:plength){
   p1max <- c(p1max, max(P1[i,]))
   p1mean <- c(p1mean, mean(P1[i,]))
   p1sum <- c(p1sum, sum(P1[i,]))
}

Could anyone please tell me how to modify and run the script in order to use all 4 CPUs?

Answer 1

Use foreach and doMC. The detailed explanation can be found here. Your script will change very little, the line

for(i in 1:plength){

should be changed to

foreach(i=1:plength) %dopar% { 

The prerequisites for any multitasking script using these packages are

library(foreach)
library(doMC)
registerDoMC()

Note of caution. According to the documentation you cannot use this in GUI.

As for your problem, do you really need multitasking? Your data.frame takes about 1.2GB of RAM, so it should fit into your memory. So you can simply use apply:

p1smry <- apply(P1,1,summary)

The result will be a matrix with summaries of each row.

You can also use function mclapply which is in the package multicore. Then your script might look like this:

loopfun <- function(i) {
     summary(P1[i,])
}

res <- mclapply(1:nrow(P1),loopfun)

This will return the list, where i-th element will be the summary of i-th row. You can convert it to matrix using sapply

mres <- sapply(res,function(x)x)

Answer 2

You've already got an answer as to how to use more than one core, but the real problem is with the way you have written your loops. Never extend your result vector/object at each iteration of a loop. If you do this, you force R to copy your result vector/object and extend it which all takes time. Instead, preallocate enough storage space before you start the loop and fill in as you go along. Here is an example:

set.seed(1)
p1 <- matrix(rnorm(10000), ncol=100)
system.time({
p1max <- p1mean <- p1sum <- numeric(length = 100)
for(i in seq_along(p1max)){
   p1max[i] <- max(p1[i,])
   p1mean[i] <- mean(p1[i,])
   p1sum[i ]<- sum(p1[i,])
}
})

   user  system elapsed 
  0.005   0.000   0.005

Or you can do these things via apply():

system.time({
p1max2 <- apply(p1, 1, max)
p1mean2 <- apply(p1, 1, mean)
p1sum2 <- apply(p1, 1, sum)
})
   user  system elapsed 
  0.007   0.000   0.006 

But note that this is no faster than doing the loop properly and sometimes slower.

However, always be on the lookout for vectorised code. You can do row sums and means using rowSums() and rowMeans() which are quicker than either the loop or apply versions:

system.time({
p1max3 <- apply(p1, 1, max)
p1mean3 <- rowMeans(p1)
p1sum3 <- rowSums(p1)
})

   user  system elapsed 
  0.001   0.000   0.002 

If I were a betting man, I would have money on the third approach I mention beating foreach() or the other multi-core options in a speed test on your matrix because they would have to speed things up considerably to justify the overhead incurred in setting up the separate processes that are farmed out the the different CPU cores.

Update: Following the comment from @shabbychef is it faster to do the sums once and reuse in the computation of the mean?

system.time({
    p1max4 <- apply(p1, 1, max)
    p1sum4 <- rowSums(p1)
    p1mean4 <- p1sum4 / ncol(p1)
    })

   user  system elapsed 
  0.002   0.000   0.002

Not in this test run, but this is far from exhaustive...

Answer 3

Have a look at the snow and snowfall packages. Plenty of examples with those...

If you want to speed up that specific code rather than learning about R and parallelism you should do this

P1 = matrix(rnorm(1000), ncol=10, nrow=10
apply(P1, 1, max)
apply(P1, 1, mean)
apply(P1, 1, sum)