r - 当方程相互依赖时，替代 for 循环进行快速计算

Question

我正在使用 for-loop进行逐步计算，其中几个方程相互依赖。由于这种依赖性，我找不到在数据帧内进行计算的解决方案。我的主要动机是在 Time 时加快计算速度。在下面的 reprex 中，向量非常大。
您能否提出以下建议的替代方案for-loop基于计算，最好在 R 中的数据帧内？我唯一能想到的就是使用 for-loop在 Rcpp .
可重现的例子

last_time <- 10
STEP = 1

Time <- seq(from = 0, to = last_time, by = STEP)

## empty vectors 
eq1 <- vector(mode = "double", length = length(Time)) 
eq2 <- vector(mode = "double", length = length(Time)) 
eq <- vector(mode = "double", length = length(Time)) 
eq3 <- vector(mode = "double", length = length(Time)) 
eq4 <- vector(mode = "double", length = length(Time)) 

## adding the first values
eq1[1] <- 25
eq2[1] <- 25
eq[1] <- 25 
eq3[1] <- 100
eq4[1] <- 2

for (t in 2:length(Time)) { 
  ## eq1
  eq1[t] <- eq[t-1] + (2.5 *  STEP * (1 - (eq[t-1])/25))
  ## eq2
  eq2[t] <- (-2 * STEP) + ((-2^2) * (STEP^2)) - (2 * eq3[t-1]) - (eq[t-1] * STEP)
  ## min.
  eq[t] <- min(eq1[t], eq2[t] )
  ## eq3
  eq3[t] <- (eq[t] - eq[t-1])/(STEP)
  ## eq4
  eq4[t] <- eq4[t-1] + (eq[t-1] * STEP) + (0.5 * eq3[t-1] * (STEP)^2)
}

输出:

my_data <- data.frame(Time, eq1, eq2, eq, eq3, eq4) 
my_data
#>    Time        eq1        eq2         eq        eq3        eq4
#> 1     0   25.00000   25.00000   25.00000 -256.00000     2.0000
#> 2     1   25.00000 -231.00000 -231.00000   25.60000  -101.0000
#> 3     2 -205.40000  225.00000 -205.40000   23.04000  -319.2000
#> 4     3 -182.36000  199.40000 -182.36000   20.73600  -513.0800
#> 5     4 -161.62400  176.36000 -161.62400   18.66240  -685.0720
#> 6     5 -142.96160  155.62400 -142.96160   16.79616  -837.3648
#> 7     6 -126.16544  136.96160 -126.16544   15.11654  -971.9283
#> 8     7 -111.04890  120.16544 -111.04890   13.60489 -1090.5355
#> 9     8  -97.44401  105.04890  -97.44401   12.24440 -1194.7819
#> 10    9  -85.19961   91.44401  -85.19961   11.01996 -1286.1037
#> 11   10  -74.17965   79.19961  -74.17965    0.00000 -1365.7934

创建于 2021-02-28 由 reprex package (v1.0.0)

Answer 1

你可以定义一个递归函数。 loop is faster than recursion尽管。

g <- function(m, STEP, time, x=2) {
  if (time == 0) m
  else {
    ## eq1
    m[x, 2] <- m[x - 1, 1] + 2.5*STEP*(1 - (m[x - 1, 1])/25)
    ## eq2
    m[x, 3] <- -2*STEP + -2^2*STEP^2 - 2*m[x - 1, 4] - m[x - 1, 1]*STEP
    ## min.
    m[x, 1] <- min(m[x, 2], m[x, 3])
    ## eq3
    m[x - 1, 4] <- (m[x, 1] - m[x - 1, 1])/STEP
    ## eq4
    m[x, 5] <- m[x - 1, 5] + m[x - 1, 1]*STEP + 0.5*m[x - 1, 4]*STEP^2
    g(m, STEP, time - 1, x + 1)
  }
}

用法

last_time <- 10; STEP <- 1

First <- c(eq0=25, eq1=25, eq2=25, eq3=100, eq4=2)
m <- matrix(0, last_time + 1, length(First), dimnames=list(NULL, names(First)))
m[1, ] <- First
    
g(m, STEP, last_time)
#              eq0        eq1        eq2        eq3        eq4
#  [1,]   25.00000   25.00000   25.00000 -256.00000     2.0000
#  [2,] -231.00000   25.00000 -231.00000   25.60000  -101.0000
#  [3,] -205.40000 -205.40000  225.00000   23.04000  -319.2000
#  [4,] -182.36000 -182.36000  199.40000   20.73600  -513.0800
#  [5,] -161.62400 -161.62400  176.36000   18.66240  -685.0720
#  [6,] -142.96160 -142.96160  155.62400   16.79616  -837.3648
#  [7,] -126.16544 -126.16544  136.96160   15.11654  -971.9283
#  [8,] -111.04890 -111.04890  120.16544   13.60489 -1090.5355
#  [9,]  -97.44401  -97.44401  105.04890   12.24440 -1194.7819
# [10,]  -85.19961  -85.19961   91.44401   11.01996 -1286.1037
# [11,]  -74.17965  -74.17965   79.19961    0.00000 -1365.7934