Эффективно стандартизировать большую матрицу в R с помощью foreach/doParallel?

Мне нужно стандартизировать (вычесть среднее значение и разделить на стандартное отклонение) столбцы нескольких больших матриц в R (примерно 300 000 строк по 10 000–20 000 столбцов).

Процесс был очень медленным, поэтому я попытался ускорить его с помощью foreach/doParallel, но, похоже, это не помогло.

Несколько примеров кода с меньшей матрицей:

library('doParallel')
n.cores <- 7
clust <- makeCluster(n.cores,type='FORK')
registerDoParallel(cl = clust)

rows <- 300000
cols <- 1500
big_matrix <- matrix(runif (rows*cols),nrow=rows,ncol=cols,dimnames=list(paste0('r',1:rows),paste0('c',1:cols)))

# Attempt 1: No parallelization. Took 52.3 seconds.
system.time( stand1 <- scale(big_matrix,center=T,scale=T) )

# Attempt 2: Each column separately with foreach. Took 168.9 seconds.
system.time(
  stand2 <- foreach(mchunk=big_matrix,.combine='cbind') %dopar% { #
              (mchunk - mean(mchunk))/sd(mchunk)
            }
)

# Attempt 3: In chunks with foreach. Took 52.7 seconds.
stand_chunks <- function(mat,n.cores) {
  chunk_indices <- split(1:ncol(mat),cut(1:ncol(mat),n.cores)) # A list of column indices, splitting the matrix into n.cores chunks
  stand_mat <- foreach(ci=chunk_indices,.combine='cbind') %dopar% {
    scale(mat[,ci],center=T,scale=T)
  }
  return(stand_mat)
}
system.time( stand3 <- stand_chunks(big_matrix,n.cores) )

# Attempt 4: In chunks with foreach, without copying whole matrix to each worker. 23 seconds.
stand_chunks_v2 <- function(mat,n.cores) {
  chunk_indices <- split(1:ncol(mat),cut(1:ncol(mat),n.cores)) # A list of column indices, splitting the matrix into n.cores chunks
  mchunks <- lapply(chunk_indices,function(ci) { return(mat[,ci]) })
  stand_mat <- foreach(mc=mchunks,.combine='cbind') %dopar% {
    scale(mc,center=T,scale=T)
  }
  return(stand_mat)
}
system.time( stand4 <- stand_chunks_v2(big_matrix,n.cores) )


parallel::stopCluster(cl=clust)

# Attempt 5: Do matrix algebra and let BLAS handle it. 20.5 seconds.
RhpcBLASctl::blas_set_num_threads(7)
stand_matalg <- function(x) {
  means <- rep(1, nrow(x)) %*% t(colMeans(x))
  sds <- rep(1, nrow(x)) %*% t(apply(x,2,sd))
  return((x - means)/sds)
}
system.time( stand5 <- stand_matalg(big_matrix) )
summary(c(abs(stand1 - stand5))) # Still looks correct, plus or minus a bit of floating-point weirdness

Я могу понять, почему попытка 2 была медленной: затраты на настройку рабочих процессов только для того, чтобы каждый из них выполнял одну быструю небольшую операцию, перевешивают преимущества распараллеливания.

Попытка 3 была очень неэффективна с использованием памяти - проблема описана здесь как «предоставление работникам большего количества данных, чем им нужно» - если я наблюдаю за потоками с помощью команды «top», я вижу, что использование памяти каждым потоком увеличивается как вся большая матрица копируется для каждого работника. К тому времени, когда я увеличиваю матрицу до 5000 столбцов, это копирование становится настолько расточительным, что ему не хватает памяти на вычислительном узле HPC с 80 ГБ ОЗУ. Это также было не быстрее, чем без распараллеливания, возможно, из-за накладных расходов на копирование большой матрицы для каждого исполнителя.

Попытка 4, наконец, оказалась быстрее. Это определенно более эффективно с точки зрения использования памяти, чем попытка 2 — она дублирует большую матрицу при составлении списка mchunks, но, по крайней мере, дублирование лучше, чем копирование ее на 7 рабочих процессов. Но, вероятно, еще можно улучшить это дальше.

Попытка 5 (адаптированная из этого ответа на старый вопрос) на данный момент является самой быстрой, но она по-прежнему требует создания двух больших матриц «средних» и «sds», которые соответствуют размерам исходной большой матрицы.

Какой более эффективный способ стандартизировать эти большие матрицы (с помощью foreach или другого метода)?