[Vegan-commits] r2213 - in pkg/vegan: R man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Thu Jun 7 10:23:17 CEST 2012


Author: jarioksa
Date: 2012-06-07 10:23:17 +0200 (Thu, 07 Jun 2012)
New Revision: 2213

Modified:
   pkg/vegan/R/simper.R
   pkg/vegan/man/simper.Rd
Log:
Squashed commit of parallel processing in simper:

commit 1d9f4d271633bec61f361292539cc189912217a8
Merge: 68018d8 e6403d4
Author: jarioksa <jari.oksanen at oulu.fi>
Date:   Thu Jun 7 01:13:25 2012 -0700

    Merge pull request #12 from EDiLD/dev-simper

    dev-simper

commit e6403d4a9abf49a7061d142e8906e276c75bf197
Author: edisz <szoe8822 at uni-landau.de>
Date:   Sun Jun 3 18:37:24 2012 +0200

    flatten loop...

commit 51ab649b6562635006657f11eda7aa7d09309e91
Author: edisz <szoe8822 at uni-landau.de>
Date:   Sun Jun 3 18:02:27 2012 +0200

    updated documentation

commit cce42b35c957f63cf61dacf28502ca35b2387890
Author: edisz <szoe8822 at uni-landau.de>
Date:   Sun Jun 3 17:46:13 2012 +0200

    added warning (empty rows)

commit e849723aeceb81ce61e893c1dad2e134dee84163
Author: edisz <szoe8822 at uni-landau.de>
Date:   Sun Jun 3 17:42:46 2012 +0200

    parallel processing of permutations

commit ef0809513af9a8aaee1db0a0cc4709ccbd8038f9
Merge: 18d2ec3 6412186
Author: edisz <szoe8822 at uni-landau.de>
Date:   Sun Jun 3 03:33:00 2012 +0200

    Merge branch 'dev-simper' of github.com:EDiLD/vegan into dev-simper

commit 18d2ec3247fd4d69c3249003b405cf4e48cd0b98
Author: edisz <szoe8822 at uni-landau.de>
Date:   Thu Apr 26 13:14:15 2012 +0200

    work on parallel processing
    updated documentation
    small speed improvement

commit 64121860b390b4ca3c9278c41ba52054571ac4c4
Author: edisz <szoe8822 at uni-landau.de>
Date:   Thu Apr 26 13:14:15 2012 +0200

    work on parallel processing
    updated documentation
    small speed improvement

Modified: pkg/vegan/R/simper.R
===================================================================
--- pkg/vegan/R/simper.R	2012-06-03 09:11:50 UTC (rev 2212)
+++ pkg/vegan/R/simper.R	2012-06-07 08:23:17 UTC (rev 2213)
@@ -1,6 +1,22 @@
 `simper` <-
-    function(comm, group, permutations = 0, trace = FALSE,  ...)
+    function(comm, group, permutations = 0, trace = FALSE,  
+             parallel = getOption("mc.cores"), ...)
 {
+    if (any(rowSums(comm, na.rm = TRUE) == 0)) 
+        warning("you have empty rows: results may be meaningless")
+    pfun <- function(x, comm, comp, i, contrp) {
+        groupp <- group[perm[x,]]
+        ga <- comm[groupp == comp[i, 1], ] 
+        gb <- comm[groupp == comp[i, 2], ]
+        for(j in 1:n.b) {
+            for(k in 1:n.a) {
+                mdp <- abs(ga[k, ] - gb[j, ])
+                mep <- ga[k, ] + gb[j, ]
+                contrp[(j-1)*n.a+k, ] <- mdp / sum(mep)  
+            }
+        }
+        colMeans(contrp)
+    }
     comm <- as.matrix(comm)
     comp <- t(combn(unique(as.character(group)), 2))
     outlist <- NULL
@@ -21,6 +37,17 @@
     nperm <- nrow(perm)
     if (nperm > 0)
         perm.contr <- matrix(nrow=P, ncol=nperm)
+    ## Parallel processing ?
+    if (is.null(parallel) && getRversion() >= "2.15.0")
+        parallel <- get("default", envir = parallel:::.reg)
+    if (is.null(parallel) || getRversion() < "2.14.0")
+        parallel <- 1
+    hasClus <- inherits(parallel, "cluster")
+    isParal <- (hasClus || parallel > 1) && require(parallel)
+    isMulticore <- .Platform$OS.type == "unix" && !hasClus
+    if (isParal && !isMulticore && !hasClus) {
+        parallel <- makeCluster(parallel)
+    }
     for (i in 1:nrow(comp)) {
         group.a <- comm[group == comp[i, 1], ]
         group.b <- comm[group == comp[i, 2], ]
@@ -36,24 +63,26 @@
         }
         average <- colMeans(contr)
         
+        ## Apply permutations
         if(nperm > 0){
             if (trace)
                 cat("Permuting", paste(comp[i,1], comp[i,2], sep = "_"), "\n")
             contrp <- matrix(ncol = P, nrow = n.a * n.b)
-            for(p in 1:nperm){
-                groupp <- group[perm[p,]]
-                ga <- comm[groupp == comp[i, 1], ] 
-                gb <- comm[groupp == comp[i, 2], ]
-                for(j in 1:n.b) {
-                    for(k in 1:n.a) {
-                        mdp <- abs(ga[k, ] - gb[j, ])
-                        mep <- ga[k, ] + gb[j, ]
-                        contrp[(j-1)*n.a+k, ] <- mdp / sum(mep)  
-                    }
-                }
-                perm.contr[ ,p] <- colMeans(contrp)
+
+            if (isParal) {
+                if (isMulticore){
+                    perm.contr <- mclapply(1:nperm, function(d) 
+                        pfun(d, comm, comp, i, contrp), mc.cores = parallel)
+                    perm.contr <- do.call(cbind, perm.contr)
+                } else {
+                    perm.contr <- parSapply(parallel, 1:nperm, function(d) 
+                        pfun(d, comm, comp, i, contrp))
+                }  
+            } else {
+                perm.contr <- sapply(1:nperm, function(d) 
+                    pfun(d, comm, comp, i, contrp))
             }
-        p <- (apply(apply(perm.contr, 2, function(x) x >= average), 1, sum) + 1) / (nperm + 1)
+            p <- (rowSums(apply(perm.contr, 2, function(x) x >= average)) + 1) / (nperm + 1)
         } 
         else {
           p <- NULL
@@ -71,6 +100,9 @@
                     avb = avb, ord = ord, cusum = cusum, p = p)
         outlist[[paste(comp[i,1], "_", comp[i,2], sep = "")]] <- out
     }
+    ## Close socket cluster if created here
+    if (isParal && !isMulticore && !hasClus)
+        stopCluster(parallel)
     attr(outlist, "permutations") <- nperm
     class(outlist) <- "simper"
     outlist
@@ -95,7 +127,9 @@
     function(object, ordered = TRUE, digits = max(3, getOption("digits") - 3), ...)
 {
     if (ordered) {
-        out <- lapply(object, function(z) data.frame(contr = z$average, sd = z$sd, ratio = z$ratio, av.a = z$ava, av.b = z$avb)[z$ord, ])
+        out <- lapply(object, function(z) 
+            data.frame(contr = z$average, sd = z$sd, ratio = z$ratio, 
+                       av.a = z$ava, av.b = z$avb)[z$ord, ])
         cusum <- lapply(object, function(z) z$cusum)
         for(i in 1:length(out)) {
             out[[i]]$cumsum <- cusum[[i]]
@@ -105,7 +139,9 @@
         } 
     } 
     else {
-        out <- lapply(object, function(z) data.frame(cbind(contr = z$average, sd = z$sd, 'contr/sd' = z$ratio, ava = z$ava, avb = z$avb, p = z$p)))
+        out <- lapply(object, function(z) 
+            data.frame(cbind(contr = z$average, sd = z$sd, 'contr/sd' = z$ratio, 
+                             ava = z$ava, avb = z$avb, p = z$p)))
     }
     attr(out, "digits") <- digits
     attr(out, "permutations") <- attr(object, "permutations")

Modified: pkg/vegan/man/simper.Rd
===================================================================
--- pkg/vegan/man/simper.Rd	2012-06-03 09:11:50 UTC (rev 2212)
+++ pkg/vegan/man/simper.Rd	2012-06-07 08:23:17 UTC (rev 2213)
@@ -10,7 +10,8 @@
 }
 
 \usage{
-simper(comm, group, permutations = 0, trace = FALSE, ...)
+simper(comm, group, permutations = 0, trace = FALSE, 
+    parallel = getOption("mc.cores"), ...)
 \method{summary}{simper}(object, ordered = TRUE,
     digits = max(3,getOption("digits") - 3), ...)
 }
@@ -29,6 +30,11 @@
   \item{ordered}{Logical; Should the species be ordered by their
     average contribution?}
   \item{digits}{Number of digits in output.}
+  \item{parallel}{Number of parallel processes or a predefined socket
+    cluster.  With \code{parallel = 1} uses ordinary, non-parallel
+    processing. The parallel processing is done with \pkg{parallel}
+    package which is available only for \R 2.14.0 and later. See 
+    \code{\link{vegandocs}} \code{decision-vegan} for details.}
   \item{...}{Parameters passed to other functions. In \code{simper} the
     extra parameters are passed to \code{\link[permute]{shuffleSet}} if
     permutations are used.}



More information about the Vegan-commits mailing list