Nested for loops in R using foreach function and doParallel library

Question

I am trying to calculate the cosine similarity between columns in a matrix. I am able to get it to work using standard for loops, but when I try to make it run in parallel to make the code run faster it does not give me the same answer. The problem is that I am unable to get the same answer using the foreach loop approach. I suspect that I am not using the correct syntax, because I have had single foreach loops work. I have tried to make the second loop a regular for loop and I used the %:% parameter with the foreach loop, but then the function doesn't even run.

Please see my attached code below. Thanks in advance for any help.

## Function that calculates cosine similarity using paralel functions.

#for calculating parallel processing
library(doParallel)

## Set up cluster on 8 cores

cl = makeCluster(8)

registerDoParallel(cl)

#create an example data
x=array(data=sample(1000*100), dim=c(1000, 100))

## Cosine similarity function using sequential for loops

cosine_seq =function (x) {

  co = array(0, c(ncol(x), ncol(x)))

  for (i in 2:ncol(x)) {
    for (j in 1:(i - 1)) {

      co[i, j] = crossprod(x[, i], x[, j])/sqrt(crossprod(x[, i]) * crossprod(x[, j]))
    }
  }

  co = co + t(co)

  diag(co) = 1

  return(as.matrix(co))

}

## Cosine similarity function using parallel for loops

cosine_par =function (x) {

  co = array(0, c(ncol(x), ncol(x)))

  foreach (i=2:ncol(x)) %dopar% {

    for (j in 1:(i - 1)) {

      co[i, j] = crossprod(x[, i], x[, j])/sqrt(crossprod(x[, i]) * crossprod(x[, j]))
    }
  }

  co = co + t(co)

  diag(co) = 1

  return(as.matrix(co))

}

## Calculate cosine similarity

tm_seq=system.time(
{

  x_cosine_seq=cosine_seq(x)

})

tm_par=system.time(
{

  x_cosine_par=cosine_par(x)

})

## Test equality of cosine similarity functions

all.equal(x_cosine_seq, x_cosine_par)

#stop cluster
stopCluster(cl)

J_F · Accepted Answer · 2017-01-31 14:22:52Z

4

The correct parallelization of nested loop uses %:% (read here).

library(foreach)
library(doParallel)
registerDoParallel(detectCores())    
cosine_par1 <- function (x) {  
  co <- foreach(i=1:ncol(x)) %:%
    foreach (j=1:ncol(x)) %dopar% {    
      co = crossprod(x[, i], x[, j])/sqrt(crossprod(x[, i]) * crossprod(x[, j]))
    }
  matrix(unlist(co), ncol=ncol(x))
}

I recommend you to write it in Rcpp, rather than running it in parallel, because foreach(i=2:n, .combine=cbind) will not always bind the columns in the right order. Also, in the above code I removed only the lower triangular condition, but the running time is considerably slower than the unparallelized code time.

set.seed(186)
x=array(data=sample(1000*100), dim=c(1000, 100))
cseq <- cosine_seq(x)
cpar <- cosine_par1(x)
 all.equal(cpar, cseq)
#[1] TRUE
head(cpar[,1])
#[1] 1.0000000 0.7537411 0.7420011 0.7496145 0.7551984 0.7602620
head(cseq[,1])
#[1] 1.0000000 0.7537411 0.7420011 0.7496145 0.7551984 0.7602620

Addendum: For this specific problem, (semi-) vectorization of cosine_seq is possible; cosine_vec is about 40-50 times faster than cosine_seq.

cosine_vec <- function(x){
  crossprod(x) / sqrt(tcrossprod(apply(x, 2, crossprod)))
}
all.equal(cosine_vec(x), cosine_seq(x))
#[1] TRUE
library(microbenchmark)
microbenchmark(cosine_vec(x), cosine_seq(x), times=20L, unit="relative")
#Unit: relative
#          expr      min       lq     mean   median       uq      max neval
# cosine_vec(x)  1.00000  1.00000  1.00000  1.00000  1.00000  1.00000    20
# cosine_seq(x) 55.81694 52.80404 50.36549 52.17623 49.56412 42.94437    20

edited Jan 31, 2017 at 14:22

J_F

10.4k2 gold badges35 silver badges59 bronze badges

answered Jan 19, 2015 at 18:48

Khashaa

7,3732 gold badges23 silver badges37 bronze badges

Sign up to request clarification or add additional context in comments.

3 Comments

Kevin Over a year ago

Thanks for the help, but I am still unable to validate the parallel representation using the "all.equal" R function. When I use the %:% syntax the cosine_par code actually only gives the identity matrix while the cosine_seq gives the correct matrix using the cosine similarity. Also I am looking at Rcpp functions, but it is not easily clear to me how to do that. Perhaps you could attach the code that does this and shows that the 2 outputs are equivalent? Thanks.

Kevin Over a year ago

Thanks the cosine_vec function is much faster! I was trying to make it even faster by using mclapply instead of apply, since I will be working with much larger matrices than this 1000x100 example. The cosine_par1 function you created does have the same output and it is much slower. However, it does not need to go over all i,j combinations since the output is symmetric (that is why I used the lower triangular condition). Won't the foreach function work with my original i,j for loop condition statements, and then add the lower triangular condition?

Yorgos Over a year ago

If there you need to insert some commands between the two "foreach", how do you do it? I get an error like "%:%" was passed an illegal right operand

Yu Yang · Accepted Answer · 2020-05-01 19:47:13Z

To do nested loop in foreach and use parallel implementation, there are two ways.

%:% + %dopar%
%dopar% + %do%

Note that for (1), it actually creates one foreach object, and you cannot add anything in between. Otherwise, you will have an error message: "%:%" was passed an illegal right operand.

And for (2), you can insert whatever you want in between. But do remember to add foreach to the .package argument in the outer loop, since the inner foreach uses foreach package.

The following is a neat way to solve the cosine matrix problem. Note that to illustrate (2), I added one extra line, and please remember to remove it for the cosine matrix calculation.

testfunc <- function (x) {
  cl<-makeCluster(4)
  registerDoParallel(cl)
  co <- foreach(i=1:ncol(x), .combine = 'rbind', .packages = c('foreach', 'stats')) %dopar% {
    k <- rnorm(3)
    foreach (j=1:ncol(x), .combine = 'c') %do% {
      crossprod(x[, i], x[, j])/sqrt(crossprod(x[, i]) * crossprod(x[, j])) + k - k
    }
  }
  stopCluster(cl)
  co
}
x <- array(data=sample(20*10), dim=c(20, 10))
testfunc(x)

Collectives™ on Stack Overflow

Nested for loops in R using foreach function and doParallel library

2 Answers 2

3 Comments

Comments

Your Answer

Linked

Hot Network Questions

Collectives™ on Stack Overflow

2 Answers 2

3 Comments

Comments

Your Answer

Sign up or log in

Post as a guest

Linked

Related