Dear list members,
Any help on this efficiency issue would be greatly appreciated.
I would like to find the most efficient way to run a non-vectorized function
(here: fisher exact test p-value) iteratively using 4 matrices with identical
dimensions. And as a result I aim for an array with identical dimensions
containing the corresponding p-values. Please consider some code using a
trivial example with 3x4 arrays below. Eventually I would like to run code on
2e3 x 7e6 arrays, for which someone suggested Amazon EC2 already...
Q1: would you agree that fisher.test() is not vectorizable? e.g. fisher.test(
matrix(c(Ax,Ay,Bx,By),ncol=2) ) does not work
Q2: direct use of Ax, Ay, Bx, By as input instead of a (list) transform for the
input would seem beneficial for speed
Q3: parallelization of the iterative process seems to make sense.
Q4: a progress bar seems to save peace of mind having no clue of the runtime.
Q5: avoidance of an output transform to get array from vector
Q6: for Q2/3/4/5 plyr seems to be ideal (e.g. maply)
Please also find some solutions below.
solution 1: using mapply
solution 2: using lapply
solution 3: using mclapply
attempt 4: stuck on plyr implementation
--Philip
### CODE START ###
Ax <- matrix(c(2,3,5,6,
3,7,8,9,
8,2,1,3), ncol = 4)
Ay <- matrix(c(9,8,5,7,
4,9,9,9,
8,7,5,4), ncol = 4)
Bx <- matrix(c(1,5,9,8,
4,7,8,9,
2,3,2,1), ncol = 4)
By <- matrix(c(5,5,9,9,
9,8,8,9,
5,5,3,2), ncol = 4)
### solution 1 using mapply
# proper answer, no input transform, output transform, no parallelization, no
progress update
sol1 <- function() {
res1 <- mapply(
function(i,j,k,l) { fisher.test( matrix(c(i,j,k,l), ncol=2),
conf.int=F)$p.value },
i=Ax, j=Ay, k=Bx, l=By,
SIMPLIFY=TRUE)
ans1 <- matrix(res1,ncol=4)
return(ans1)
}
s1 <- sol1()
### solution 2 using lapply
# proper answer, input transform, output transform, no parallelization, no
progress update
sol2 <- function() {
tmp.list <- as.data.frame(rbind(as.numeric(Ax), as.numeric(Ay),
as.numeric(Bx), as.numeric(By)))
# determine fisher.test p-values as list
res2 <- lapply(tmp.list,
function(x) { fisher.test( matrix(unlist(x), ncol=2), conf.int=F)$p.value
})
ans2 <- matrix(unlist(res2),ncol=4)
return(ans2)
}
s2 <- sol2()
### solution 3 using mclapply
# proper answer using input transform, output transform, parallelization, no
progress update
library(multicore)
sol3 <- function() {
tmp.list <- as.data.frame(rbind(as.numeric(Ax), as.numeric(Ay),
as.numeric(Bx), as.numeric(By)))
# determine fisher.test p-values as list
res3 <- mclapply(tmp.list,
function(x) { fisher.test( matrix(unlist(x), ncol=2), conf.int=F)$p.value },
mc.cores=4)
ans3 <- matrix(unlist(res3),ncol=4)
}
s3 <- sol3()
### solution 4 using plyr::maply
# difficulty finding equivalent code
# benefit could be: no input transform, no output transform, parallelization,
and progress update
library(plyr)
library(abind)
library(doMC)
registerDoMC(cores=4)
sol4 <- function() {
ans4 <- maply(
#.data = abind(i=Ax,j=Ay,k=Bx,l=By,along=0),
#.data = abind(Ax,Ay,Bx,By,along=3),
#.data = data.frame(i=Ax, j=Ay, k=Bx, l=By),
#.data = cbind(i=as.vector(Ax), j=as.vector(Ay), k=as.vector(Bx),
l=as.vector(By)),
#.data = list(i=Ax, j=Ay, k=Bx, l=By),
.data = abind(i=Ax, j=Ay, k=Bx, l=By, along=3),
.fun = function(i,j,k,l) { fisher.test( matrix(c(i,j,k,l), ncol=2),
conf.int=F)$p.value },
.progress = "text",
.parralel = TRUE
)
return(ans4)
}
all.equal(s1,s2) # TRUE
all.equal(s1,s3) # TRUE
library(microbenchmark)
microbenchmark(sol1, times=1000)
microbenchmark(sol2, times=1000)
microbenchmark(sol3, times=1000)
microbenchmark(sol4, times=1000)
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
