Re: [問題] 加速模擬時的矩陣相乘

看板R_Language作者 (攸藍)時間10年前 (2014/07/07 19:54), 10年前編輯推噓0(000)
留言0則, 0人參與, 最新討論串2/2 (看更多)
改寫之後可以快上200倍左右 (主要差異是Reduce跟apply的運算方式不同,可以自行google。) 如果可以利用C++,直接寫迴圈可以更快,這裡我就不獻醜了。 我用比較少的值做計算,如果是您的size的話,應該加速不止200倍 另外,其中 Fsim = F[m, ] + eps 我覺得您需要確定一下是否與您所想的相同 M = 3 N = 500 T = 100 nSim = 300 # 稍微減少一下 B = matrix(rnorm(N * M), nrow = N) F = matrix(rnorm(M * T), nrow = M) sigma = abs(rnorm(M)) t1 = Sys.time() Z2 = lapply(1:M, function(m){ eps = matrix(rnorm(nSim * T, 0, sigma[m]), nrow = nSim) Fsim = t(apply(eps, 1, function(v) v + F[m,])) lapply(1:nSim, function(n)outer(B[,m], Fsim[n,])) }) result2 = Reduce(function(x, y){ lapply(1:length(x), function(i) x[[i]] + y[[i]]) }, Z2) t2 = Sys.time() t2 - t1 # Time difference of 0.2330141 secs t1 = Sys.time() Z = array(0, dim = c(M, nSim, N, T)) # preallocate for (m in 1:M) { eps = rnorm(nSim * T, 0, sigma[m]) #generate random normal dim(eps) = c(nSim, T) # 變成 10000 x T Fsim = F[m, ] + eps for (n.path in 1:nSim) { Z[m,n.path,,] = as.matrix(B[,m]) %*% Fsim[n.path,] } } result = apply(Z, c(2,3,4), sum) t2 = Sys.time() t2 - t1 # Time difference of 52.60001 secs 如果要parallel,就library(snowfall) 可以參考拙作 #1Hocf0cY (R_Language) ※ 引述《kolun (...)》之銘言: : 我有一個模型Y = BF : 其中B是N x M, F是M x T的矩陣 : 我要透過模擬F的變化10000次再組合回去得到模擬的Y : 但速度非常慢 : 我現在的作法是 : # data : M = 3 : N = 500 : T = 100 : nSim = 10000 : B = matrix(rnorm(N * M), nrow = N) : F = matrix(rnorm(M * T), nrow = M) : sigma = abs(rnorm(M)) : Z = array(0, dim = c(M, nSim, N, T)) # preallocate : for (m in 1:M) { : eps = rnorm(nSim * T, 0, sigma[m]) #generate random normal : dim(eps) = c(nSim, T) # 變成 10000 x T : Fsim = F[m, ] + eps : for (n.path in 1:nSim) { : Z[m,n.path,,] = as.matrix(B[,m]) %*% Fsim[n.path,] : } : } : result = apply(Z, c(2,3,4), sum) : 其中第一步allocation就要大概3秒 (on macbook air 2013 mid) : 產生亂數的時間好像還在預期中 : 但是裡面這個for (npath in 1:10000) {...}的時間實在太久了 : 我知道我現在的作法應該是最慢的 : 想請教版上的各位大大有什麼方法可以提升速度呢? : 謝謝各位 -- ※ 發信站: 批踢踢實業坊(ptt.cc), 來自: 140.116.152.225 ※ 文章網址: http://www.ptt.cc/bbs/R_Language/M.1404734080.A.91E.html ※ 編輯: celestialgod (218.164.196.107), 07/08/2014 00:06:39
文章代碼(AID): #1Jkeg0aU (R_Language)
討論串 (同標題文章)
文章代碼(AID): #1Jkeg0aU (R_Language)