Comparing mclapply vs foreach

Cleanup our outcome data

gmvol<-readNIfTI('../../neuroimaging/4d6mmnew.nii.gz')
mask<-readNIfTI('../../neuroimaging/vol0000.nii.gz')
ymat <- ymat_trans_mask(neurodat = gmvol, mask = mask)

Cleanup our predictor variables

xmat <- data.frame(read_excel('../../neuroimaging/VMB Variables.xlsx',skip=1))
xtx <- xmat_trans(xmat)

Comparing Permutation test between mclapply vs foreach

permtest.oneunit<-function(xmat=xtx,y=ymat,seed=1234){
  set.seed(seed)
  n=nrow(xtx)
  indx=sample(n,n,replace=FALSE)
  bhat<- y %*% xmat[indx, 2:3] #####-----why only second and third column (groups), test the beta coefficient of the groups
  return(bhat)
}

run_pertest<-function(x.mat=xtx,y=ymat){
  
  vindx=rep(1:101,each=floor(nrow(ymat)/100))[1:nrow(ymat)]
  
  bhat = y %*% x.mat[,2:3]
  
  rr1=lapply(1:101, function(k){
    perms<-mclapply(1:100,function(seednum){permtest.oneunit(xmat=xtx,seed=seednum,y=y[vindx ==k,])},mc.cores=8)
    aa1=do.call(rbind,lapply(perms,function(x)x[,1]))
    ps1=apply(rbind(bhat[vindx==k,1],aa1),2, function(x)sum( x[-1]> abs(x[1]) |  x[-1]< (-1*abs(x[1])) ))/100
    
    aa2=do.call(rbind,lapply(perms,function(x)x[,2]))
    ps2=apply(rbind(bhat[vindx==k,2],aa2),2, function(x)sum( x[-1]> abs(x[1]) |  x[-1]< (-1*abs(x[1])) ))/100
 
    return(cbind(unlist(ps1),unlist(ps2)))})
  
  return(do.call(rbind,rr1))
}

Runtime using mclapply

system.time(run1 <- run_pertest(x.mat=xtx,y=ymat))

##    user  system elapsed 
##   9.817   6.008  28.969

Runtime using foreach

cl <- makeCluster(8)
registerDoParallel(cl)
system.time(run2 <- run_permtest_multicat(x.mat = xtx, y = ymat, columns = 2:3, split = 101, num_perms = 100))

##    user  system elapsed 
##   5.248   0.437  25.033

Verifying same results

print(cbind(head(run1), head(run2)))

##                result.1 result.2
## [1,] 0.49 0.88     0.49     0.88
## [2,] 0.60 0.76     0.60     0.76
## [3,] 0.44 0.88     0.44     0.88
## [4,] 0.51 0.79     0.51     0.79
## [5,] 0.59 0.73     0.59     0.73
## [6,] 0.72 0.66     0.72     0.66

As we can see, by utilizing a foreach structure, built into the run_permtest_culticat function from this package, rather than that of one that relies on mclapply, we have cut around 20% of runtime for permutation test calculations. Performance may vary depending on operating system (this was done on macOS).