Комментарии 12
Первое правило оптимизаций — не оптимизировать то, что и так хорошо работает.
Начать нужно с примитивных функций, например с sum.of.squares
и лишнего apply
.
Вот пример того, что я получаю. Тестирую apply
по разным измерениям, затем purrr::map_dbl
(самый быстрый) и future::map_dbl
(ооооочень медленно).
library(tidyverse)
library(microbenchmark)
#> Warning: package 'microbenchmark' was built under R version 4.0.3
library(furrr)
#> Warning: package 'furrr' was built under R version 4.0.3
#> Loading required package: future
#> Warning: package 'future' was built under R version 4.0.3
a <- matrix(rnorm(500000, mean=0, sd=2), 100000, 50)
sum_of_squares <- function(n) {
n_sq <- n ^ 2
sum(n_sq)
}
# Verifying summing over different dimensions
abs(
sum(apply(a, 1, sum_of_squares)) -
sum(apply(a, 2, sum_of_squares)))
#> [1] 0
abs(
sum(map_dbl(1:dim(a)[2], ~sum_of_squares(a[,.x]))) -
sum(map_dbl(1:dim(a)[1], ~sum_of_squares(a[.x, ]))))
#> [1] 0
plan(cluster, workers = 4L)
# Verifying cluster setup, should be 4 unique pids
future_map_int(1:20, ~Sys.getpid()) %>% unique
#> [1] 14240 2112 5356 14848
abs(
sum(future_map_dbl(1:dim(a)[2], ~sum_of_squares(a[,.x]))) -
sum(future_map_dbl(1:dim(a)[1], ~sum_of_squares(a[.x, ]))))
#> [1] 0
microbenchmark(
apply_dim_1 = sum(apply(a, 1, sum_of_squares)),
apply_dim_2 = sum(apply(a, 2, sum_of_squares)),
map_dim_1 = sum(map_dbl(1:dim(a)[1], ~sum_of_squares(a[.x, ]))),
map_dim_2 = sum(map_dbl(1:dim(a)[2], ~sum_of_squares(a[,.x]))),
future_dim_1 = sum(future_map_dbl(1:dim(a)[1], ~sum_of_squares(a[.x, ]))),
future_dim_2 = sum(future_map_dbl(1:dim(a)[2], ~sum_of_squares(a[,.x]))),
times = 100L)
#> Unit: milliseconds
#> expr min lq mean median uq max neval
#> apply_dim_1 404.4141 460.1293 520.03837 510.4109 563.45995 802.6375 100
#> apply_dim_2 130.5241 150.7894 192.64263 172.8856 234.17170 553.5947 100
#> map_dim_1 445.9119 476.0938 528.86206 516.6661 549.29145 913.7160 100
#> map_dim_2 39.1026 53.3201 60.41777 60.3456 68.07625 86.4652 100
#> future_dim_1 789.9283 877.1860 954.96807 936.8079 1012.73945 1566.2939 100
#> future_dim_2 601.0447 660.5598 720.77519 692.2518 752.43105 1506.9873 100
Created on 2020-11-09 by the reprex package (v0.3.0)
Поправьте, если где-то ошибка. При этом все тесты запускались на подручном лэптопе на Core-i5 с двумя физическими ядрами.
Самое очевидное — нет абсолютно никакого смысла считать вдоль dim = 1
, даже если с точки зрения memory layout это верно. map_dim_2
выполняет суммирование всех квадратов вдоль dim = 2
за 60 миллисекунд, вызывая кастомную функцию только 50 раз, что почти в 9 раз быстрее чем вдоль dim = 1
, что приводит к вызову той же функции 100_000 раз.
Опять же, если я нигде не допустил ошибки, использование multi-process parallelization в данном случае (на жалких 5 миллионах точек) не оправдано. С другой стороны, если каждый из 50 вызовов суммирующей функции занимал бы, например, 10 секунд, то распараллеливание дало бы ускорение, на вскидку, порядка ~2 раз на моем CPU, и наверное близко к ~10 на CPU автора (в зависимости от оверхеда).
Предложенная идея с подсчётом по столбцам в данном случае совсем не подходит: так можно поступать, только если мы хотим просуммировать все значения в матрице, но нам нужны отдельные результаты для каждой из строк. Так что никакие способы с перегруппировкой вычислений не походят. Грубо говоря, три велосипеда – это три руля, три рамы и шесть колёс. Но если одному ребёнку дать в подарок три руля, второму – три рамы, а третьему – шесть колёс, это будет для каждого из них совсем не то же, что получить в подарок целый велосипед. Хотя суммы и равны. Допустим, если строки – это показания с набора датчиков на каждый момент времени, то вычисление суммы квадратов показаний с первого датчика за всю историю наблюдений никак не поможет нам получить сумму квадратов показаний всех датчиков в тот или иной момент времени.
А, ну я значит упустил момент про абсолютную наобходимость суммирования по строкам. Если величины настолько разные, почему тогда матрица а не data.frame
? Есть какие-то ограничения? Работать с прямоугольными данными проще, если они — data.frame
.
В любом случае, я лишь хотел отметить, что просто так добавить future
— не панацея.
Например, в вашем случае у вас прирост примерно в 4 раза, хотя логических ядер у вас 12. Т.е., предполагая что availableWorkers()
у вас тоже 12, можно еще что-то оптимизировать прежде чем распараллеливать. Возможно, меньшее количество процессов даст больший прирост производительности на данном конкретрном демо-примере.
Отдельный вопрос, полезный для изучения — как future_apply*
распределяет работу?
Например, функция может отгружать новые задачи по мере завершения предыдущих, по одной. Альтернативно — сразу разделить весь датасет и назначить каждому процессу свой блок, после чего ожидать завершения всех участвующих процессов. В первом случае большое количество легковесных операций (как в данном примере) затормозит работу, и, возможно, потребует дополнительных оптимизаций.
Уменьшение числа используемых ядер (по умолчанию используются все ядра) выглядит не самой хорошей идеей. Непонятно, каким образом это может повысить производительность. В моём случае наилучшие результаты были с параметрами по умолчанию. А насчёт возможностей тонкой настройки параметров параллельного запуска согласен — там вполне может быть ещё какой-то запас пространства для оптимизаций. Единственное, это уже скорее всего будет специфично конкретно для данной аппаратной конфигурации и решаемой задачи. Но попробовать можно.
Ну вот с таким детальным разъяснением становится понятно, зачем демо-пример именно в таком виде. Если вы применили такой подход на реальных данных, можете поделиться результатом? Хотя бы порядок ускорения на N ядер (как у вас тут, 3.8 раз на 12 логических ядер)?
Но в любом случае, есть какой-то предел, выше которого подняться в конкретном случае не получается. При работе с числами ускорение на 12 ядрах у меня в лучшем случае обычно приближается к 4. Судя по тому, что при преобразовании матрицы в таблицу данных ускорение стремится уже не к 4, а к 8 (при существенно общей потере производительности), 4 – это не какой-то абсолютный предел, и на данных других типов ускорение может быть и выше.
На моём ПК в данном примере самая быстрая реализация оказалась в пакете matrixStats,
хотя есть ошибка на пределе машинной точности (от -9.095e-13 до 9.095e-13).
library(matrixStats)
start_time <- Sys.time()
a2 <- a^2
bRowSums <- rowSums2(a2)
timediff <- difftime(Sys.time(), start_time)
cat("Расчёт занял: ", timediff, units(timediff))
# Расчёт занял: 0.01400399 secs
identical(b, bRowSums)
# FALSE
summary(b - bRowSums)
# Min. 1st Qu. Median Mean 3rd Qu. Max.
# -9.095e-13 -2.842e-14 0.000e+00 8.440e-17 2.842e-14 9.095e-13
В data.table расчёт дольше, хотя при работе с десятками миллионов записей широких таблиц data.table меня много раз выручал.
library(data.table)
dt <- as.data.table(a)
start_time <- Sys.time()
dt2 <- dt^2
colNames <- names(dt)
dt2[, SUM := rowSums(.SD), .SDcols = colNames]
timediff <- difftime(Sys.time(), start_time)
cat("Расчёт занял: ", timediff, units(timediff))
# Расчёт занял: 0.05101109 secs
identical(dt2$SUM, b)
# TRUE
library(tictoc)
tic()
for(i in 1:dim(a)[1]) {
b[i] <- sum.of.squares(a[i,])
}
toc()
# 4.16 sec elapsed
library(doParallel)
cl <- makeCluster(16)
registerDoParallel(cl)
tic()
parB <- foreach(i=seq(dim(a)[1]), .inorder = FALSE, .combine = rbind) %dopar% {
sum.of.squares(a[i,])
}
toc()
# 23.62 sec elapsed
При этом в более пригодном формате хранения больших таблиц (в data.table) пример считается в разы быстрее, если нам действительно нужно что-то построчно сделать:
library(data.table)
library(matrixStats)
dt <- as.data.table(a)
dt2 <- dt^2
colNames <- names(dt)
tic()
dt2[, SUM := rowSums(.SD), .SDcols = colNames]
toc()
# 0.03 sec elapsed
identical(dt2$SUM, b)
# TRUE
Реализация в пакете future.apply у меня в разы медленнее оказалась на 3950х:
library("future.apply")
plan(multisession)
tic()
bFuture <- future_apply(a, 1, function(x) sum.of.squares(x))
toc()
# 2.61 sec elapsed
Вообще, когда доходит до вычислительных задач, то обилие пакетов в R из преимущества становится скорее проблемой. Почти никто не пользуется чистым R, все ставят пакеты. Пакетов тысячи, их функции пересекаются и дублируют друг друга, как эти функции соотносятся в плане производительности — во многих случаях просто неизвестно, некоторые вещи никто никогда не сравнивал. У каждого из пользователей — своя уникальная экосистема из десятков привычных ему пакетов, то есть фактически вместо единого языка существует огромное число частично пересекающихся между собой диалектов. Можно стать экспертом по одному диалекту, можно владеть несколькими, но практически невозможно охватить их все. Из-за этого теряется какая-либо возможность сделать путеводитель по правильным техникам решения тех или иных задач в R, а изучение возможностей языка начинает носить характер натурных экспериментов и какой-то алхимии. Иногда это даже забавно, но когда дело касается вычислений — весёлого мало. Хочется идти по правильному пути, писать сразу более-менее начисто — а официального правильного пути нет, есть тысячи пересекающихся тропинок и надо самому сравнивать их между собой. Это становится реальной проблемой, слишком уж много времени уходит на сравнения различных способов. Многие советуют переписывать тяжёлые части на каком-то более производительном языке и вызывать их из R, возможно в некотором смысле это самое правильное решение, хотя и достаточно громоздкое. Пока что для меня это открытая проблема.
Такой вариант в лоб, который дает тот же результат,
rowSums(a*a)
вычисляется на моей машине в среднем за 0.016 секунд. Среднее время расчета вашего исходного примера с циклом у меня занимает 4.43 секунды.
На мой взгляд, пример из виньетки future.apply тоже прост и показателен — средствами базового R здесь future.apply не превзойти.
x <- replicate(1e5, runif(100), simplify = F)
y <- lapply(x, quantile, probs = 1:3/4)
yy <- future_lapply(x, quantile, probs = 1:3/4)
all.equal(y, yy)
Распараллеливаем код в R за пару минут