Skip to content
Browse files

added util functions; turned off parallelization

  • Loading branch information...
1 parent 843ef1e commit f49e4b101550e07cda9597488001e1cc5ab5d7ed Michael Milham committed Oct 5, 2012
Showing with 90 additions and 31 deletions.
  1. +61 −0 simulations/lib_utils.R
  2. +29 −31 simulations/main.R
View
61 simulations/lib_utils.R
@@ -4,3 +4,64 @@ zdist <- function(mat, method="cor") {
stop("Unsupported method for zdist")
)
}
+
+vcat <- function(verbose, msg, ..., newline=TRUE) {
+ if (verbose) {
+ cat(sprintf(msg, ...))
+ if (newline) cat("\n")
+ }
+}
+
+vstop <- function(msg, ...) stop(sprintf(msg, ...))
+
+vsystem <- function(verbose, cmd, ...) {
+ vcat(verbose, cmd, ...)
+ ret <- system(sprintf(cmd, ...))
+ if (ret != 0)
+ stop("command failed")
+}
+
+set_parallel_procs <- function(nforks=1, nthreads=1, verbose=FALSE) {
+ is.installed <- function(mypkg) is.element(mypkg, installed.packages()[,1])
+
+ vcat(verbose, "Setting %i parallel forks", nforks)
+ suppressPackageStartupMessages(library("doMC"))
+ registerDoMC()
+ nprocs <- getDoParWorkers()
+ if (nforks > nprocs) {
+ vstop("# of forks %i is greater than the actual # of processors (%i)",
+ nforks, nprocs)
+ }
+ options(cores=nforks)
+
+ vcat(verbose, "Setting %i threads for matrix algebra operations",
+ nthreads)
+ #nprocs <- omp_get_max_threads()
+ if (nthreads > nprocs) {
+ vstop("# of threads %i is greater than the actual # of processors (%i)",
+ nthreads, nprocs)
+ }
+
+ if (existsFunction("setMKLthreads", where=topenv(.GlobalEnv))) {
+ vcat(verbose, "...using Intel's MKL")
+ setMKLthreads(nthreads)
+ } else if (is.installed("blasctl")) {
+ # cover all our blases
+ vcat(verbose, "...using GOTOBLAS or Other")
+ suppressPackageStartupMessages(library("blasctl"))
+ if (existsFunction("blas_set_num_threads", where=topenv(.GlobalEnv)))
+ blas_set_num_threads(nthreads)
+ omp_set_num_threads(nthreads)
+ } else {
+ vcat(verbose, "...NOT using any multithreading BLAS library")
+ }
+
+ # Not sure if these env vars matter?
+ Sys.setenv(OMP_NUM_THREADS=nthreads)
+ Sys.setenv(GOTO_NUM_THREADS=nthreads)
+ Sys.setenv(MKL_NUM_THREADS=nthreads)
+ Sys.setenv(OMP_DYNAMIC=TRUE)
+ Sys.setenv(MKL_DYNAMIC=TRUE)
+
+ invisible(TRUE)
+}
View
60 simulations/main.R
@@ -6,26 +6,20 @@ source("lib_mdmr.R")
source("lib_utils.R")
# Settings
-run_parallel <- TRUE
-ncores <- 3
-nthreads <- 2
+run_parallel <- FALSE # getting segfaults???
+nforks <- 3
+nthreads <- 4
nsubs <- 200
nnodes <- 100
nnei <- 10
effect_sizes <- c(0.02, 0.04, 0.06, 0.08, 0.1, 0.15, 0.2, 0.3)
effect_size_noise <- 0.01 # 1%
num_nodes_change <- nnodes/20
num_conns_change_per_node <- c(1:5, 6, 8, 10, 12, 15)
+verbose <- TRUE
-if (run_parallel) {
- library(doMC)
- registerDoMC(ncores)
-
- if (system("hostname", intern=T) == "rocky") {
- library(blasctl)
- omp_set_num_threads(nthreads)
- }
-}
+if (run_parallel)
+ set_parallel_procs(nforks, nthreads, verbose)
grp <- factor(rep(c("A","B"),each=nsubs/2))
@@ -36,12 +30,14 @@ grp <- factor(rep(c("A","B"),each=nsubs/2))
###
# Generate group adjacency matrix
+vcat(verbose, "Generating group adjacency matrix")
g <- watts.strogatz.game(1, nnodes, nnei, 0.05)
g <- simplify(g) # removes loops and multiple connections
group.adj <- get.adjacency(g)
# TODO: HERE I DON'T TAKE INTO ACCOUNT THAT THIS IS A UNDIRECTED GRAPH
# Add weights
+vcat(verbose, "Adding weights")
group.wt <- group.adj
conns <- which(group.adj==1)
no_conns <- which(group.adj==0)
@@ -62,17 +58,19 @@ while (TRUE) {
group.wt[no_conns][bad] <- rnorm(sum(bad), 0, 0.1)
}
+vcat(verbose, "")
###
# Subject Connectivity Matrices
###
# each subject is a random variation on the group average
+vcat(verbose, "Generating subject weighted matrices")
subjects.wt <- laply(1:nsubs, function(i) {
subject.wt <- group.wt
subject.wt[conns] <- subject.wt[conns] + runif(length(conns), min=-0.1, max=0.1)
- subject.wt
-}, .progress="text", .parallel=F)
+ as.matrix(subject.wt)
+}, .progress="text", .parallel=run_parallel)
# 1. create group difference (i.e., add effect size)
subjects_with_diffs <- laply(effect_sizes, function(es) {
@@ -105,20 +103,20 @@ subjects_with_diffs <- laply(effect_sizes, function(es) {
})
}, .progress="text", .parallel=run_parallel)
-
-
-###
-# Distance Matrices
-###
-
-distances_with_diffs <- llply(distances_with_diffs, function(Dss) {
-
- llply(Dss, function(Ds) {
- aaply(Ds) {
- as.vector(zdist(t(xs[,,ni])))
- })
- })
-
- aaply(, 3, function())
-
-})
+#
+#
+####
+## Distance Matrices
+####
+#
+#distances_with_diffs <- llply(distances_with_diffs, function(Dss) {
+#
+# llply(Dss, function(Ds) {
+# aaply(Ds) {
+# as.vector(zdist(t(xs[,,ni])))
+# })
+# })
+#
+# aaply(, 3, function())
+#
+#})

0 comments on commit f49e4b1

Please sign in to comment.
Something went wrong with that request. Please try again.