Последовательный, совокупный расчет

Мне нужно сделать расчет временных рядов, где значение, рассчитанное в каждой строке, зависит от результата, вычисленного в предыдущей строке. Я надеюсь использовать удобствоdata.table, Актуальной проблемой является гидрологическая модель - вычисление совокупного водного баланса, добавление осадков на каждом временном шаге и вычитание стока и испарения в зависимости от текущего объема воды. Набор данных включает в себя различные бассейны и сценарии (группы). Здесь я буду использовать более простую иллюстрацию проблемы.

Упрощенный пример расчета выглядит так, для каждого временного шага (строки)i:

 v[i] <- a[i] + b[i] * v[i-1]

a а такжеb являются векторами значений параметров, иv является вектором результата. Для первого ряда (i == 1) начальная стоимостьv принимается какv0 = 0.

Первая попытка

Моей первой мыслью было использоватьshift() вdata.table, Минимальный пример, включая желаемый результатv.ans, является

library(data.table)        # version 1.9.7
DT <- data.table(a = 1:4, 
                 b = 0.1,
                 v.ans = c(1, 2.1, 3.21, 4.321) )
DT
#    a   b v.ans
# 1: 1 0.1 1.000
# 2: 2 0.1 2.100
# 3: 3 0.1 3.210
# 4: 4 0.1 4.321

DT[, v := NA]   # initialize v
DT[, v := a + b * ifelse(is.na(shift(v)), 0, shift(v))][]
#    a   b v.ans v
# 1: 1 0.1 1.000 1
# 2: 2 0.1 2.100 2
# 3: 3 0.1 3.210 3
# 4: 4 0.1 4.321 4

Это не работает, потому чтоshift(v) дает копию оригинального столбцаv, сдвинут на 1 ряд. Это не зависит от назначенияv.

Я также рассмотрел построение уравнения с использованием cumsum () и cumprod (), но это тоже не сработает.

Метод грубой силы

Поэтому для удобства я прибегаю к циклу for внутри функции:

vcalc <- function(a, b, v0 = 0) {
  v <- rep(NA, length(a))      # initialize v
  for (i in 1:length(a)) {
    v[i] <- a[i] + b[i] * ifelse(i==1, v0, v[i-1])
  }
  return(v)
}

Эта накопительная функция прекрасно работает с data.table:

DT[, v := vcalc(a, b, 0)][]
#    a   b v.ans     v
# 1: 1 0.1 1.000 1.000
# 2: 2 0.1 2.100 2.100
# 3: 3 0.1 3.210 3.210
# 4: 4 0.1 4.321 4.321
identical(DT$v, DT$v.ans)
# [1] TRUE

Мой вопрос

Мой вопрос, могу ли я написать этот расчет более кратким и эффективнымdata.table Кстати, без необходимости использовать цикл for и / или определение функции? С помощьюset() возможно?

Или есть лучший подход все вместе?

Изменить: лучший цикл

Решение Дэвида Rcpp ниже вдохновило меня удалитьifelse() отfor цикл:

vcalc2 <- function(a, b, v0 = 0) {
  v <- rep(NA, length(a))
  for (i in 1:length(a)) {
    v0 <- v[i] <- a[i] + b[i] * v0
  }
  return(v)
}

vcalc2() на 60% быстрее, чемvcalc().