使用并行优化数组 R 的循环

Optimize loops for arrays R using parallel

我有一个数组 data = array[1:50,1:50,1:50] 里面的值是 -1, 1 之间的实数。

"Data" 可以视为 50x50x50 立方体。

我需要根据这个等式创建一个相关矩阵(删除所有零)=>

值 = (x+y)-|x-y|并且矩阵大小是可能组合的 2 倍 (50x50x50)*((50x50x50)-1)/2 = 7.812.437.500 这 2 倍 = 相关矩阵。

我这样做了:

假设我们有 3x3x3:

arr = array(rnorm(10), dim=c(3,3,3))

data = data.frame(array(arr))


data$voxel <- rownames(data) 

#remove zeros
data<-data[!(data[,1]==0),]

rownames(data) = data$voxel

data$voxel = NULL


#######################################################################################
#Create cluster

no_cores <- detectCores() #- 1

clus <- makeCluster(no_cores)

clusterExport(clus, list("data") , envir=environment())

clusterEvalQ(clus,
             compare_strings <- function(j,i) {
               value <- (data[i,]+data[j,])-abs(data[i,]- data[j,])   
               pair <- rbind(rownames(data)[j],rownames(data)[i],value)
               return(pair)
             })

i = 0 # start 0
kk = 1
table <- data.frame()

ptm <- proc.time()

while(kk<nrow(data)) {

  out <-NULL  
  i = i+1 # fix row
  j = c((kk+1):nrow(data)) # rows to be compared

  #Apply the declared function  
  out = matrix(unlist(parRapply(clus,expand.grid(i,j), function(x,y) compare_strings(x[1],x[2]))),ncol=3, byrow = T)

  table <- rbind(table,out)

  kk = kk +1

}

proc.time() - ptm

结果是data.frame:

v1  v2  v3
1   2   2.70430114250358
1   3   0.199941717684129
... up to 351 rows

但这需要几天的时间...

另外我想为这个相关性创建一个矩阵:

   1                         2              3...
1  1                  2.70430114250358 
2  2.70430114250358          1
3...

有没有更快的方法?

谢谢

您的代码中存在一些性能错误:

  1. 当你应该依赖向量化时你循环了。
  2. 你在循环中增长一个对象。
  3. 您并行化了循环的每次迭代,而不是并行化了外部循环。

如果你避免了第一个问题,你就可以避免所有这些问题。

显然,您想要比较行的每个组合。为此,您应该首先获取行索引的所有组合:

combs <- t(combn(1:27, 2))

然后你可以将比较函数应用到这些:

compare <- function(j,i, data) {
  as.vector((data[i,]+data[j,])-abs(data[i,]- data[j,]))
}

res <- data.frame(V1 = combs[,1], V2 = combs[,2], 
                  V3 = compare(combs[,1], combs[,2], data))

现在,如果我们想检查这是否给出与您的代码相同的结果,我们首先需要修复您的输出。通过将字符(行名)与矩阵中的数字组合,您将得到一个字符矩阵,并且最终 data.frame 的列都是字符。之后我们可以使用 type.convert 来修复它(尽管从一开始就应该避免):

table[] <- lapply(table, function(x) type.convert(as.character(x)))

现在我们可以看到结果是一样的:

all.equal(res, table)
#[1] TRUE

如果愿意,可以将结果转成稀疏矩阵:

library(Matrix)
m <- sparseMatrix(i = res$V1, j = res$V2, x = res$V3, 
                  dims = c(27, 27), symmetric = TRUE)
diag(m) <- 1