在6個嵌套for循環代碼中節省時間

Question

我正在嘗試爲我的代碼找到正確的參數，這是擴展卡爾曼濾波器的回溯。我有6個嵌套循環，每個參數一個。目前，當我對每個參數有3個可能的值時，代碼平均需要大約5分鐘才能運行，但看到隨着我增加參數的數量，我所花費的時間會增加n^6。我有點擔心。

有沒有什麼辦法可以優化代碼以節省更多時間？

PS - 只需使用任何數據文件而不是給定的Reddy.csv（1180行數據） PPS - 最後我需要找到最小值爲i，j，k，l，m，n的值MSE。

下面的代碼：

start.time <- Sys.time() 

library(invgamma) 
w = read.csv("Reddy.csv") 
q = ts(w[2]) 
num = length(q) 

f = function(x){ 
    f1 = sqrt(x) 
    return(f1) 
} 
h = function(x){ 
    h1 = x**3 
    return(h1) 
} 


ae1 = seq(24,26,0.1) 
ae2 = seq(24,26,0.1) 

be1 = seq(0.1,2,0.1) 
be2 = seq(0.1,2,0.1) 

a = seq(1,3,0.1) 
b = seq(0.1,2,0.1) 

count = 0 

MSE = matrix(nrow = length(ae1)*length(ae2)*length(be1)*length(be2)*length(a)*length(b), ncol =7) 

for (i in ae1){ 
    for (j in ae2){ 
    for (k in be1){ 
     for (l in be2){ 
     for (m in a){ 
      for (n in b){ 
      d = rep(0,num) 
      xt = rep(0,num) 
      yt = rep(0,num) 
      fx = rep(0,num) 
      hx = rep(0,num) 

      e = rinvgamma(num,i,k) 
      g = rinvgamma(num,j,l) 
      for(o in 2:num){ 
       fx[o] = f(xt[o-1]) 
       xt[o] = m*fx[o] + e[o-1] 
       hx[o] = h(xt[o]) 
       yt[o]= n*hx[o] +g[o] 
       d[o] = (yt[o] - q[o])**2 
      } 
      count <- count + 1 
      MSE[count,1] = mean(d) 
      MSE[count,2] = i 
      MSE[count,3] = j 
      MSE[count,4] = k 
      MSE[count,5] = l 
      MSE[count,6] = m 
      MSE[count,7] = n 
      t = rbind(mean(d),i,j,k,l,m,n) 
      print(t) 
      } 
     } 
     } 
    } 
    } 
} 

end.time <- Sys.time() 
time.taken <- end.time - start.time 
time.taken 

m = which.min(MSE[,1]) 
MSE[m,] 

Answer 1

進一步優化可以通過矢量化一些這樣的最內層循環的計算來實現，

start.time <- Sys.time() 

library(invgamma) 
w = rnorm(1180) 
q = ts(w) 
num = length(q) 

f = function(x){ 
    f1 = sqrt(x) 
    return(f1) 
} 
h = function(x){ 
    h1 = x**3 
    return(h1) 
} 


ae1 = seq(24,26,1) 
ae2 = seq(24,26,1) 

be1 = seq(0.1,2,0.7) 
be2 = seq(0.1,2,0.7) 

a = seq(1,3,1) 
b = seq(0.1,2,0.7) 

count = 0 

MSE = matrix(nrow = length(ae1)*length(ae2)*length(be1)*length(be2)*length(a)*length(b), ncol =7) 

for (i in ae1){ 
    for (j in ae2){ 
     for (k in be1){ 
      for (l in be2){ 
       for (m in a){ 
        for (n in b){ 
         d = rep(0,num) 
         xt = rep(0,num) 
         yt = rep(0,num) 
         fx = rep(0,num) 
         hx = rep(0,num) 
         e = rinvgamma(num,i,k) 
         g = rinvgamma(num,j,l) 
         for(o in 2:num){ 
          fx[o] = f(xt[o-1]) 
          xt[o] = m*fx[o] + e[o-1] 
         } 
         hx = h(xt) 
         yt = n*hx +g 
         d = (yt - q)**2 

         count <- count + 1 
         MSE[count,1] = mean(d) 
         MSE[count,2] = i 
         MSE[count,3] = j 
         MSE[count,4] = k 
         MSE[count,5] = l 
         MSE[count,6] = m 
         MSE[count,7] = n 
         ## t = rbind(mean(d),i,j,k,l,m,n) 
         ## print(t) 
        } 
       } 
      } 
     } 
    } 
} 

end.time <- Sys.time() 
time.taken <- end.time - start.time 
time.taken 

m = which.min(MSE[,1]) 
MSE[m,] 

在我的筆記本電腦，這可以使代碼幾次更快。

編輯：

如果組合是太多了，你不希望有大的矩陣創建和修改，因爲你只關心最小MSE的，你只能不停的最佳戰績像這樣的組合：

MSEopt <- Inf 

for (i in ae1){ 
    for (j in ae2){ 
     for (k in be1){ 
      for (l in be2){ 
       for (m in a){ 
        for (n in b){ 
         d = rep(0,num) 
         xt = rep(0,num) 
         yt = rep(0,num) 
         fx = rep(0,num) 
         hx = rep(0,num) 
         e = rinvgamma(num,i,k) 
         g = rinvgamma(num,j,l) 
         for(o in 2:num){ 
          fx[o] = f(xt[o-1]) 
          xt[o] = m*fx[o] + e[o-1] 
         } 
         hx = h(xt) 
         yt = n*hx +g 
         d = (yt - q)**2 

         if (mean(d) < MSEopt) { 
          ## print(MSEopt) 
          MSEopt <- mean(d) 
          best_combn <- list(i = i, j = j, k = k, l = l, m = m, n = n) 
         } 
        } 
       } 
      } 
     } 
    } 
}