Если верить стереотипам, то язык R – это что-то узкоспециализированное для статистики и машинного обучения. Второй стереотип – что код на чистом R не очень быстрый: во-первых, потому что интерпретируемый, во-вторых, потому что исполняется последовательно. Безусловно, стереотипы имеют какую-то связь с реальностью, иначе бы их не существовало, но на то они и стереотипы, что дают экстремально упрощённую картину мира, в которой теряется много деталей. В частности, сегодня я хочу поделиться удивительно простым способом, как добавить параллельность в R и кратно ускорить выполнение уже имеющегося кода без необходимости внесения в него каких-либо серьёзных изменений. Всё это делается буквально за пару минут.
Допустим, у нас есть матрица или таблица данных, содержащая некоторое количество строк и столбцов, и мы хотим для каждой из строк выполнить какое-то однотипное вычисление. Например, посчитать сумму квадратов её значений. Логично вынести расчёты в функцию и вызывать её для каждой из строк.
Исходные данные:
a <- matrix(rnorm(500000, mean=0, sd=2), 100000, 50)
Функция:
sum.of.squares <- function(n) {
n_sq <- sapply(n, function(x) x^2)
sum(n_sq)
}
Можно просто перебрать циклом строки, и к каждой из строк применить эту функцию, но это не самый рекомендуемый способ для R. Вычисления для каждой из строк будут выполняться последовательно, все расчёты будут идти на одном ядре. Такой код действительно не очень производителен. На всякий случай запишем такой вариант и замерим время его выполнения:
b <- vector()
for(i in 1:dim(a)[1]) {
b[i] <- sum.of.squares(a[i,])
}
Замеряем время выполнения:
b <- vector()
start_time <- Sys.time()
for(i in 1:dim(a)[1]) {
b[i] <- sum.of.squares(a[i,])
}
timediff <- difftime(Sys.time(), start_time)
cat("Расчёт занял: ", timediff, units(timediff))
Получаем:
Расчёт занял: 4.474074 secs
Будем использовать это время как некоторую отправную точку для сравнения с другими способами.
Выполнение заданной операции над каждой из строк некоторого объекта можно записать в векторизованной форме. Для этого в R есть специальное семейство функций apply(). Для тех кто не знаком с ними, могу порекомендовать статьи по их использованию: 1, 2. Данные функции примечательны тем, что позволяют компактно описать однотипные действия без использования циклов. Даже в нашей функции расчёта суммы квадратов я для краткости использовал одну из таких функций – sapply(), которая выполняет некоторое заданное действие для каждого из элементов поступающего на вход набора значений. В нашем случае – вычисляет квадрат каждого из значений. Если же говорить про применение данной функции к матрице исходных данных, то с использованием функции apply() это можно переписать как:
b <- apply(a, 1, function(x) sum.of.squares(x))
Действительно, очень компактная запись. Но если замерить время её исполнения, то получим примерно то же, что и у цикла:
start_time <- Sys.time()
b <- apply(a, 1, function(x) sum.of.squares(x))
timediff <- difftime(Sys.time(),start_time)
cat("Расчёт занял: ", timediff, units(timediff))
Результат:
Расчёт занял: 4.484046 secs
К слову, уже первый знак после запятой может от запуска к запуску немного меняться. Но суть ясна: это примерно тот же результат, что и с циклом.
Казалось бы, мы возвращаемся к исходному положению, что R считает всё в одном потоке, и от этого значительная часть ресурсов компьютера простаивает без дела. В некотором смысле так и есть: в данный момент функции семейства apply(), хотя потенциально их и можно распараллелить, выполняются последовательно. Но если рассмотреть ситуацию немного шире, то уже сейчас идёт работа, чтобы сделать их параллельными. Более того, уже сейчас эти функции из будущего можно загрузить и использовать вместо обычных функций семейства apply(). Помимо непосредственно функции apply() распараллеливаются также функции by(), eapply(), lapply(), Map(), .mapply(), mapply(), replicate(), sapply(), tapply(), vapply(). Со временем параллельные реализации этих функций заменят текущие, а пока что для их использования необходимо поставить специальный пакет future_apply:
install.packages("future.apply")
Буквально несколько секунд – и вот он установлен. После этого его надо подключить и указать, что расчёт пойдёт в многоядерном режиме:
library("future.apply")
plan(multiprocess)
Есть различные параметры запуска. Например, поддерживается расчёт на распределённом кластере. Текущие настройки параллельного расчёта можно посмотреть используя команду future::plan(). Чтобы пользоваться возможностями параллельного расчёта, код основной программы менять практически никак не нужно, достаточно дописать к функциям apply приставку "future_". Получим:
b <- future_apply(a, 1, function(x) sum.of.squares(x))
Запустим с замером времени:
start_time <- Sys.time()
b <- future_apply(a, 1, function(x) sum.of.squares(x))
timediff <- difftime(Sys.time(),start_time)
cat("Расчёт занял: ", timediff, units(timediff))
Получаем уже совсем другое время выполнения:
Расчёт занял: 1.283569 secs
В моём случае расчёт производился на Intel Core i7-8750H с 12 логическими ядрами. Полученное ускорение конечно не 12-кратное, но всё же довольно приличное.
В зависимости от решаемой задачи степень ускорения может существенно меняться. В некоторых случаях, например, на маленьких массивах, накладные расходы от распараллеливания могут превышать выигрыш от параллельного исполнения, и вычисления могут даже замедляться, так что везде без разбора бездумно использовать параллельные версии функций точно не стоит. Например, не стоит распараллеливать таким способом расчёт квадратов значений внутри вызываемой функции, используя future_sapply, это приведёт к сильному замедлению. Но для обработки массивов с большим числом строк подобное распараллеливание почти всегда даёт достаточно ощутимое ускорение. В данном случае – примерно четырёхкратное. Если производить аналогичные вычисления не для матрицы, а для таблицы данных, например, предварительно преобразовав в неё исходную матрицу (a <- data.frame(a)), то расчёт в целом в несколько раз замедлится, но разница между последовательным и параллельным выполнением будет уже примерно в 8 раз. Довольно ощутимая разница.
Ну вот, собственно, и всё. Способ достаточно простой. Для меня, когда я про него узнал, это была просто находка. Справедливо ли утверждение, что нынешний R не поддерживает параллельные вычисления? Зависит от точки зрения на этот вопрос, от строгости его постановки. Но в некотором смысле можно считать, что всё-таки поддерживает.
BkmzSpb
Первое правило оптимизаций — не оптимизировать то, что и так хорошо работает.
Начать нужно с примитивных функций, например с
sum.of.squares
и лишнегоapply
.Вот пример того, что я получаю. Тестирую
apply
по разным измерениям, затемpurrr::map_dbl
(самый быстрый) иfuture::map_dbl
(ооооочень медленно).Created on 2020-11-09 by the reprex package (v0.3.0)
A_Degteryov Автор
Задача, на которой всё демонстрируется – это не реальный случай, а специально упрощённый пример, придуманный просто для удобства изложения. Обычно на вход приходят не случайные значения с нормальным распределением, а внутри функции сидит какое-то более сложное вычисление, чем просто сумма квадратов. Цель была не в том, чтобы показать, как оптимизировать конкретно этот пример, а чтобы продемонстрировать способ ускорения одинаковых операций над строками за счёт распараллеливания вычислений. При этом хотелось показать всё на максимально понятном примере, чтобы не отвлекаться на объяснения каких-то специфических вычислений.
Предложенная идея с подсчётом по столбцам в данном случае совсем не подходит: так можно поступать, только если мы хотим просуммировать все значения в матрице, но нам нужны отдельные результаты для каждой из строк. Так что никакие способы с перегруппировкой вычислений не походят. Грубо говоря, три велосипеда – это три руля, три рамы и шесть колёс. Но если одному ребёнку дать в подарок три руля, второму – три рамы, а третьему – шесть колёс, это будет для каждого из них совсем не то же, что получить в подарок целый велосипед. Хотя суммы и равны. Допустим, если строки – это показания с набора датчиков на каждый момент времени, то вычисление суммы квадратов показаний с первого датчика за всю историю наблюдений никак не поможет нам получить сумму квадратов показаний всех датчиков в тот или иной момент времени.
BkmzSpb
А, ну я значит упустил момент про абсолютную наобходимость суммирования по строкам. Если величины настолько разные, почему тогда матрица а не
data.frame
? Есть какие-то ограничения? Работать с прямоугольными данными проще, если они —data.frame
.В любом случае, я лишь хотел отметить, что просто так добавить
future
— не панацея.Например, в вашем случае у вас прирост примерно в 4 раза, хотя логических ядер у вас 12. Т.е., предполагая что
availableWorkers()
у вас тоже 12, можно еще что-то оптимизировать прежде чем распараллеливать. Возможно, меньшее количество процессов даст больший прирост производительности на данном конкретрном демо-примере.Отдельный вопрос, полезный для изучения — как
future_apply*
распределяет работу?Например, функция может отгружать новые задачи по мере завершения предыдущих, по одной. Альтернативно — сразу разделить весь датасет и назначить каждому процессу свой блок, после чего ожидать завершения всех участвующих процессов. В первом случае большое количество легковесных операций (как в данном примере) затормозит работу, и, возможно, потребует дополнительных оптимизаций.
A_Degteryov Автор
Основное преимущество таблиц данных (data.frame) — возможность хранить в разных столбцах значения различающихся типов (числовые значения, категориальные переменные, текст и т.п.). По всей видимости, эта возможность не бесплатная, поэтому если она не нужна (например, если мы обрабатываем показания с набора датчиков и это заведомо будут только числа), то использование матрицы может дать некоторый прирост производительности. Опять же, я не занимался специальной оптимизацией конкретно этого примера, он исключительно демонстрационный, но даже в нём переход с матрицы на таблицу данных даёт замедление расчёта большее, чем можно отыграть распараллеливанием.
Уменьшение числа используемых ядер (по умолчанию используются все ядра) выглядит не самой хорошей идеей. Непонятно, каким образом это может повысить производительность. В моём случае наилучшие результаты были с параметрами по умолчанию. А насчёт возможностей тонкой настройки параметров параллельного запуска согласен — там вполне может быть ещё какой-то запас пространства для оптимизаций. Единственное, это уже скорее всего будет специфично конкретно для данной аппаратной конфигурации и решаемой задачи. Но попробовать можно.
BkmzSpb
Ну вот с таким детальным разъяснением становится понятно, зачем демо-пример именно в таком виде. Если вы применили такой подход на реальных данных, можете поделиться результатом? Хотя бы порядок ускорения на N ядер (как у вас тут, 3.8 раз на 12 логических ядер)?
A_Degteryov Автор
Степень ускорения в первую очередь зависит от количества строк входных данных и «тяжести» обрабатывающей процедуры, куда в качестве составляющей в некотором смысле входит и количество элементов строки. Замедление происходит либо при слишком малом числе строк (допустим, если будет не 500000 строк, а 5000 – прирост будет уже менее чем двукратный), либо при слишком простой процедуре, требующей мало вычислений. Со случаями, где сильно влияла бы структура данных или порядок значений, я не сталкивался, так что условно можно считать, что значительного влияния нет и случайные значения в этом плане мало отличимы от реальных.
Но в любом случае, есть какой-то предел, выше которого подняться в конкретном случае не получается. При работе с числами ускорение на 12 ядрах у меня в лучшем случае обычно приближается к 4. Судя по тому, что при преобразовании матрицы в таблицу данных ускорение стремится уже не к 4, а к 8 (при существенно общей потере производительности), 4 – это не какой-то абсолютный предел, и на данных других типов ускорение может быть и выше.
Vladic82
И точно, даже в данном маленьком примере матрица отработала быстрее «быстрых» таблиц:
На моём ПК в данном примере самая быстрая реализация оказалась в пакете matrixStats,
хотя есть ошибка на пределе машинной точности (от -9.095e-13 до 9.095e-13).
В data.table расчёт дольше, хотя при работе с десятками миллионов записей широких таблиц data.table меня много раз выручал.
Vladic82
Согласен, что пример очень неудачный, но благодарю за примеры в пакете, о котором не знал. Сам привык выполнять параллельные расчёты в R в doParallel, который на данном примере вообще не успевал ядра нагрузить, но расчёт через %dopar% выполнился аж за 24 секунды, а при монопотоке как в публикации в 4 секунды:
При этом в более пригодном формате хранения больших таблиц (в data.table) пример считается в разы быстрее, если нам действительно нужно что-то построчно сделать:
Реализация в пакете future.apply у меня в разы медленнее оказалась на 3950х:
A_Degteryov Автор
Спасибо, интересно взглянуть на вопрос с другой стороны. Трюк с поэлементной обработкой значений пригоден только конкретно для этого учебного примера, и в большинстве реальных случаев будет неприменим. Попробовал переложить исходные данные в data.table, но само по себе это это дало только небольшое замедление. Возможно выигрыш будет на более крупных массивах, при случае проверю.
Вообще, когда доходит до вычислительных задач, то обилие пакетов в R из преимущества становится скорее проблемой. Почти никто не пользуется чистым R, все ставят пакеты. Пакетов тысячи, их функции пересекаются и дублируют друг друга, как эти функции соотносятся в плане производительности — во многих случаях просто неизвестно, некоторые вещи никто никогда не сравнивал. У каждого из пользователей — своя уникальная экосистема из десятков привычных ему пакетов, то есть фактически вместо единого языка существует огромное число частично пересекающихся между собой диалектов. Можно стать экспертом по одному диалекту, можно владеть несколькими, но практически невозможно охватить их все. Из-за этого теряется какая-либо возможность сделать путеводитель по правильным техникам решения тех или иных задач в R, а изучение возможностей языка начинает носить характер натурных экспериментов и какой-то алхимии. Иногда это даже забавно, но когда дело касается вычислений — весёлого мало. Хочется идти по правильному пути, писать сразу более-менее начисто — а официального правильного пути нет, есть тысячи пересекающихся тропинок и надо самому сравнивать их между собой. Это становится реальной проблемой, слишком уж много времени уходит на сравнения различных способов. Многие советуют переписывать тяжёлые части на каком-то более производительном языке и вызывать их из R, возможно в некотором смысле это самое правильное решение, хотя и достаточно громоздкое. Пока что для меня это открытая проблема.
jzha
У вас выбран не самый удачный демонстрационный пример.
вычисляется на моей машине в среднем за 0.016 секунд.Такой вариант в лоб, который дает тот же результат,
Среднее время расчета вашего исходного примера с циклом у меня занимает 4.43 секунды.
A_Degteryov Автор
Я взял максимально упрощённый пример, чтобы не нужно было устраивать лирическое отступление и рассказывать, что в нём вычисляется. При этом хотелось обойтись без псевдокода, а привести какое-то вычисление, которое действительно можно запустить и выполнить. Поскольку это только демонстрационный пример, цели максимально затюнить конкретно это вычисление не ставилось, я даже не придумывал ему защиту от решения другими способами. Допустим, завтра нам понадобится посчитать для каждой строки не сумму квадратов всех значений, а сумму квадратов разниц между текущим значением и средним значением по ряду. Послезавтра — медианное значение дельт между текущим значением и всеми остальными. Тут трюки с поэлементной обработкой уже не подойдут, а предложенный в статье способ довольно универсален, ему не очень важно, какое вычисление работает внутри
jzha
Да-да ваша цель ясна, полностью поддерживаю вас. Просто пример вызывал некоторое недоумение. Он как раз демонстрирует сильную сторону базового R.
На мой взгляд, пример из виньетки future.apply тоже прост и показателен — средствами базового R здесь future.apply не превзойти.