Ускорить цикл работы в R
У меня большая проблема с производительностью в R. Я написал функцию, которая перебирает data.frame объект. Это просто добавляет новый столбец к data.frame и что-то накапливает. (простая операция). data.frame имеет примерно 850 тыс. строк. Мой ПК все еще работает (около 10 часов), и я понятия не имею о времени выполнения.
dayloop2 <- function(temp){
for (i in 1:nrow(temp)){
temp[i,10] <- i
if (i > 1) {
if ((temp[i,6] == temp[i-1,6]) & (temp[i,3] == temp[i-1,3])) {
temp[i,10] <- temp[i,9] + temp[i-1,10]
} else {
temp[i,10] <- temp[i,9]
}
} else {
temp[i,10] <- temp[i,9]
}
}
names(temp)[names(temp) == "V10"] <- "Kumm."
return(temp)
}
Есть идеи, как ускорить эту операцию?
9 ответов
Самая большая проблема и корень неэффективности - индексирование data.frame, я имею в виду все эти строки, где вы используете temp[,],
Старайтесь избегать этого как можно больше. Я взял твою функцию, поменяй индексацию и вот version_A
dayloop2_A <- function(temp){
res <- numeric(nrow(temp))
for (i in 1:nrow(temp)){
res[i] <- i
if (i > 1) {
if ((temp[i,6] == temp[i-1,6]) & (temp[i,3] == temp[i-1,3])) {
res[i] <- temp[i,9] + res[i-1]
} else {
res[i] <- temp[i,9]
}
} else {
res[i] <- temp[i,9]
}
}
temp$`Kumm.` <- res
return(temp)
}
Как вы можете видеть, я создаю вектор res которые собирают результаты. В конце я добавляю это data.frame и мне не нужно связываться с именами. Так как же лучше?
Я запускаю каждую функцию для data.frame с nrow от 1000 до 10000 на 1000 и измерять время с system.time
X <- as.data.frame(matrix(sample(1:10, n*9, TRUE), n, 9))
system.time(dayloop2(X))
Результат
Вы можете видеть, что ваша версия экспоненциально зависит от nrow(X), Модифицированная версия имеет линейную зависимость и простую lm Модель предсказывает, что для 850 000 строк вычисление занимает 6 минут 10 секунд.
Сила векторизации
Как Шейн и Калимо заявляют в своих ответах, векторизация является ключом к повышению производительности. Из вашего кода вы можете выйти за пределы цикла:
- кондиционирование
- инициализация результатов (которые
temp[i,9])
Это приводит к этому коду
dayloop2_B <- function(temp){
cond <- c(FALSE, (temp[-nrow(temp),6] == temp[-1,6]) & (temp[-nrow(temp),3] == temp[-1,3]))
res <- temp[,9]
for (i in 1:nrow(temp)) {
if (cond[i]) res[i] <- temp[i,9] + res[i-1]
}
temp$`Kumm.` <- res
return(temp)
}
Сравните результат для этой функции, на этот раз для nrow от 10000 до 100000000.
Тюнинг настроенного
Еще один твик - изменение в цикле индексации temp[i,9] в res[i] (которые в точности повторяются в итерации i-го цикла). Это опять разница между индексацией вектора и индексацией data.frame,
Второе: когда вы смотрите на цикл, вы видите, что нет необходимости циклически повторять все i, но только для тех, которые соответствуют условиям.
Итак, поехали
dayloop2_D <- function(temp){
cond <- c(FALSE, (temp[-nrow(temp),6] == temp[-1,6]) & (temp[-nrow(temp),3] == temp[-1,3]))
res <- temp[,9]
for (i in (1:nrow(temp))[cond]) {
res[i] <- res[i] + res[i-1]
}
temp$`Kumm.` <- res
return(temp)
}
Производительность, которую вы получаете, сильно зависит от структуры данных. Именно - на процент TRUE значения в состоянии. Для моих смоделированных данных требуется время вычисления на 850 000 строк ниже одной секунды.
Если вы хотите, вы можете пойти дальше, я вижу, по крайней мере, две вещи, которые можно сделать:
- Напиши
Cкод, чтобы сделать условную сумму если вы знаете, что в ваших данных максимальная последовательность не велика, вы можете изменить цикл на векторизованное время, что-то вроде
while (any(cond)) { indx <- c(FALSE, cond[-1] & !cond[-n]) res[indx] <- res[indx] + res[which(indx)-1] cond[indx] <- FALSE }
Код, используемый для моделирования и рисунков, доступен на GitHub.
Общие стратегии для ускорения кода R
Во-первых, выясните, где медленная часть на самом деле. Нет необходимости оптимизировать код, который не работает медленно. Для небольших объемов кода, просто продумав его, можно работать. Если это не помогло, RProf и подобные инструменты профилирования могут быть полезны.
Как только вы обнаружите узкое место, подумайте о более эффективных алгоритмах для выполнения того, что вы хотите. Расчеты должны выполняться только один раз, если это возможно, поэтому:
- Сохраняйте результаты и получайте к ним доступ, а не повторно пересчитывайте
- Возьмите не зависящие от цикла вычисления из циклов
- Избегайте необязательных вычислений (например , не используйте регулярные выражения с фиксированным поиском)
Использование более эффективных функций может привести к умеренному или большому увеличению скорости. Например, paste0 дает небольшой прирост эффективности, но .colSums() и его родственники приносят более заметные выгоды. mean особенно медленно
Тогда вы можете избежать некоторых наиболее распространенных проблем:
cbindзамедлит вас очень быстро.- Инициализируйте ваши структуры данных, затем заполните их, вместо того, чтобы каждый раз расширять их.
- Даже с предварительным распределением вы можете переключиться на подход с передачей по ссылке, а не с передачей по значению, но это может не стоить хлопот.
- Взгляните на R Inferno, чтобы увидеть больше подводных камней, которых следует избегать.
Попробуйте улучшить векторизацию, которая часто, но не всегда, помогает. В связи с этим по своей сути векторизованные команды, такие как ifelse, diff и тому подобное обеспечит большее улучшение, чем apply семейство команд (которые обеспечивают практически полное повышение скорости по сравнению с хорошо написанным циклом).
Вы также можете попытаться предоставить больше информации для функций R. Например, используйте vapply скорее, чем sapply и укажите colClasses при чтении в текстовых данных. Прирост скорости будет изменяться в зависимости от того, сколько угадывания вы устраняете.
Далее рассмотрим оптимизированные пакеты: data.table Пакет может значительно увеличить скорость, где это возможно, при манипулировании данными и при чтении больших объемов данных (fread).
Затем попробуйте увеличить скорость за счет более эффективных способов вызова R:
- Скомпилируйте ваш R-скрипт. Или используйте
Raа такжеjitсогласованные пакеты для своевременной компиляции (Дирк имеет пример в этой презентации). - Убедитесь, что вы используете оптимизированный BLAS. Они обеспечивают повсеместный прирост скорости. Честно говоря, это позор, что R не использует автоматически самую эффективную библиотеку при установке. Надеюсь, что Revolution R внесет свой вклад в общее сообщество.
- Рэдфорд Нил провел ряд оптимизаций, некоторые из которых были перенесены в R Core, а многие другие были разветвлены в pqR.
И наконец, если все вышеперечисленное по-прежнему не позволяет вам работать так быстро, как вам нужно, вам может потребоваться перейти на более быстрый язык для медленного фрагмента кода. Сочетание Rcpp а также inline здесь делает замену только самой медленной части алгоритма кодом C++ особенно легко. Вот, например, моя первая попытка сделать это, и она поражает даже высоко оптимизированные R-решения.
Если после всего этого у вас все еще есть проблемы, вам просто нужно больше вычислительной мощности. Посмотрите на распараллеливание ( http://cran.r-project.org/web/views/HighPerformanceComputing.html) или даже решения на основе GPU (gpu-tools).
Ссылки на другие руководства
Если вы используете for циклы, вы, скорее всего, кодируете R, как если бы это был C или Java или что-то еще. Правильно векторизованный R-код очень быстр.
Возьмем, к примеру, эти два простых фрагмента кода для генерации списка из 10000 целых чисел в последовательности:
Первый пример кода - это то, как можно было бы закодировать цикл, используя традиционную парадигму кодирования. Это займет 28 секунд, чтобы завершить
system.time({
a <- NULL
for(i in 1:1e5)a[i] <- i
})
user system elapsed
28.36 0.07 28.61
Вы можете получить улучшение почти в 100 раз, просто выполнив предварительное выделение памяти:
system.time({
a <- rep(1, 1e5)
for(i in 1:1e5)a[i] <- i
})
user system elapsed
0.30 0.00 0.29
Но используя базовую векторную операцию R с помощью оператора двоеточия : эта операция практически мгновенная:
system.time(a <- 1:1e5)
user system elapsed
0 0 0
Это можно сделать намного быстрее, пропустив циклы с помощью индексов или вложенных ifelse() заявления.
idx <- 1:nrow(temp)
temp[,10] <- idx
idx1 <- c(FALSE, (temp[-nrow(temp),6] == temp[-1,6]) & (temp[-nrow(temp),3] == temp[-1,3]))
temp[idx1,10] <- temp[idx1,9] + temp[which(idx1)-1,10]
temp[!idx1,10] <- temp[!idx1,9]
temp[1,10] <- temp[1,9]
names(temp)[names(temp) == "V10"] <- "Kumm."
Мне не нравится переписывать код... Также, конечно, ifelse и lapply - лучшие варианты, но иногда это трудно сделать подобающим.
Часто я использую data.frames, как можно использовать списки, такие как df$var[i]
Вот вымышленный пример:
nrow=function(x){ ##required as I use nrow at times.
if(class(x)=='list') {
length(x[[names(x)[1]]])
}else{
base::nrow(x)
}
}
system.time({
d=data.frame(seq=1:10000,r=rnorm(10000))
d$foo=d$r
d$seq=1:5
mark=NA
for(i in 1:nrow(d)){
if(d$seq[i]==1) mark=d$r[i]
d$foo[i]=mark
}
})
system.time({
d=data.frame(seq=1:10000,r=rnorm(10000))
d$foo=d$r
d$seq=1:5
d=as.list(d) #become a list
mark=NA
for(i in 1:nrow(d)){
if(d$seq[i]==1) mark=d$r[i]
d$foo[i]=mark
}
d=as.data.frame(d) #revert back to data.frame
})
версия data.frame:
user system elapsed
0.53 0.00 0.53
версия списка:
user system elapsed
0.04 0.00 0.03
В 17 раз быстрее использовать список векторов, чем в data.frame.
Любые комментарии о том, почему внутренне data.frames так медленно в этом отношении? Казалось бы, они работают как списки...
Для еще более быстрого кода сделайте это class(d)='list' вместо d=as.list(d) а также class(d)='data.frame'
system.time({
d=data.frame(seq=1:10000,r=rnorm(10000))
d$foo=d$r
d$seq=1:5
class(d)='list'
mark=NA
for(i in 1:nrow(d)){
if(d$seq[i]==1) mark=d$r[i]
d$foo[i]=mark
}
class(d)='data.frame'
})
head(d)
Как сказал Ари в конце своего ответа, Rcpp а также inline пакеты позволяют невероятно легко делать вещи быстро. В качестве примера попробуйте это inline код (предупреждение: не проверено):
body <- 'Rcpp::NumericMatrix nm(temp);
int nrtemp = Rccp::as<int>(nrt);
for (int i = 0; i < nrtemp; ++i) {
temp(i, 9) = i
if (i > 1) {
if ((temp(i, 5) == temp(i - 1, 5) && temp(i, 2) == temp(i - 1, 2) {
temp(i, 9) = temp(i, 8) + temp(i - 1, 9)
} else {
temp(i, 9) = temp(i, 8)
}
} else {
temp(i, 9) = temp(i, 8)
}
return Rcpp::wrap(nm);
'
settings <- getPlugin("Rcpp")
# settings$env$PKG_CXXFLAGS <- paste("-I", getwd(), sep="") if you want to inc files in wd
dayloop <- cxxfunction(signature(nrt="numeric", temp="numeric"), body-body,
plugin="Rcpp", settings=settings, cppargs="-I/usr/include")
dayloop2 <- function(temp) {
# extract a numeric matrix from temp, put it in tmp
nc <- ncol(temp)
nm <- dayloop(nc, temp)
names(temp)[names(temp) == "V10"] <- "Kumm."
return(temp)
}
Там аналогичная процедура для #includeвещи, где вы просто передаете параметр
inc <- '#include <header.h>
к функции, как include=inc, Что действительно круто в этом, так это то, что он выполняет всю компоновку и компиляцию за вас, поэтому прототипирование действительно быстрое.
Отказ от ответственности: я не совсем уверен, что класс tmp должен быть числовым, а не числовой матрицей или чем-то еще. Но я в основном уверен.
Изменить: если вам все еще нужно больше скорости, OpenMP является средством параллелизации, подходящим для C++, Я не пробовал использовать его из inline, но это должно работать. Идея была бы, в случае n ядра, есть итерация цикла k быть выполненным k % n, Подходящее введение можно найти в Matloff's Art of R Programming, доступной здесь, в главе 16, " Использование C".
Ответы здесь великолепны. Один незначительный аспект, который не был затронут, заключается в том, что вопрос гласит: "Мой компьютер все еще работает (около 10 часов), и я понятия не имею о времени выполнения". Я всегда вставляю следующий код в циклы при разработке, чтобы понять, как изменения влияют на скорость, а также для отслеживания того, сколько времени потребуется для завершения.
dayloop2 <- function(temp){
for (i in 1:nrow(temp)){
cat(round(i/nrow(temp)*100,2),"% \r") # prints the percentage complete in realtime.
# do stuff
}
return(blah)
}
Работает с lapply также.
dayloop2 <- function(temp){
temp <- lapply(1:nrow(temp), function(i) {
cat(round(i/nrow(temp)*100,2),"% \r")
#do stuff
})
return(temp)
}
Если функция в цикле довольно быстрая, но число циклов велико, рассмотрите возможность периодической печати, так как печать на самой консоли имеет накладные расходы. например
dayloop2 <- function(temp){
for (i in 1:nrow(temp)){
if(i %% 100 == 0) cat(round(i/nrow(temp)*100,2),"% \r") # prints every 100 times through the loop
# do stuff
}
return(temp)
}
В R вы часто можете ускорить обработку цикла, используя apply семейные функции (в вашем случае это, вероятно, будет replicate). Посмотрите на plyr пакет, который предоставляет индикаторы выполнения.
Другой вариант - вообще избежать циклов и заменить их векторной арифметикой. Я не совсем уверен, что вы делаете, но вы, вероятно, можете применить свою функцию ко всем строкам одновременно:
temp[1:nrow(temp), 10] <- temp[1:nrow(temp), 9] + temp[0:(nrow(temp)-1), 10]
Это будет намного быстрее, и тогда вы можете отфильтровать строки с вашим условием:
cond.i <- (temp[i, 6] == temp[i-1, 6]) & (temp[i, 3] == temp[i-1, 3])
temp[cond.i, 10] <- temp[cond.i, 9]
Векторизованная арифметика требует больше времени и размышлений о проблеме, но иногда вы можете сэкономить несколько порядков времени выполнения.
Взгляните на
accumulate() функция от
{purrr}:
dayloop_accumulate <- function(temp) {
temp %>%
as_tibble() %>%
mutate(cond = c(FALSE, (V6 == lag(V6) & V3 == lag(V3))[-1])) %>%
mutate(V10 = V9 %>%
purrr::accumulate2(.y = cond[-1], .f = function(.i_1, .i, .y) {
if(.y) {
.i_1 + .i
} else {
.i
}
}) %>% unlist()) %>%
select(-cond)
}
Обработка с data.table это приемлемый вариант:
n <- 1000000
df <- as.data.frame(matrix(sample(1:10, n*9, TRUE), n, 9))
colnames(df) <- paste("col", 1:9, sep = "")
library(data.table)
dayloop2.dt <- function(df) {
dt <- data.table(df)
dt[, Kumm. := {
res <- .I;
ifelse (res > 1,
ifelse ((col6 == shift(col6, fill = 0)) & (col3 == shift(col3, fill = 0)) ,
res <- col9 + shift(res)
, # else
res <- col9
)
, # else
res <- col9
)
}
,]
res <- data.frame(dt)
return (res)
}
res <- dayloop2.dt(df)
m <- microbenchmark(dayloop2.dt(df), times = 100)
#Unit: milliseconds
# expr min lq mean median uq max neval
#dayloop2.dt(df) 436.4467 441.02076 578.7126 503.9874 575.9534 966.1042 10
Если вы игнорируете возможные выгоды от фильтрации условий, это очень быстро. Очевидно, что если вы можете сделать расчет на подмножестве данных, это помогает.