PTT數位生活區 / R_Language

[問題] 迴圈的計算(小複雜)

看板R_Language作者 (empireisme)時間4年前 (2019/11/25 12:34), 4年前編輯推噓1(1推 0噓 20→)
留言21則, 4人參與, 5年前最新討論串1/2 (看更多)
我本身使用R大概一兩年 目前想要寫出一個小function可以計算出以下的公式 https://imgur.com/supH7mE
d 是要計算出來的向量 長度是m A 是一個矩陣維度是n*m，裡面的元素都是正的 alpha是A矩陣的元素 p 是一個向量長度為m 其為機率向量 裡面的元素相加等於1 以下是我的程式碼，雖然檢查過很多次，但不確定有沒有算錯 或是有沒有辦法用矩陣的方式去計算之，因為這樣算太慢了 set.seed(1) m=10 n=7 A <- matrix( rexp(m*n),ncol=m ,nrow=n ) p <- c( 0.05,0.1,0.3,0.05,0.05,0.05,0.07,0.03,0.2,0.1) sum(p ) length(p ) den <- rep(0,n) densum <- rep(0,n) d<- rep(0,m) for(j in 1:m){ for(i in 1:n){ den[i]<- sum(A[i,]*p) densum[i]<- A[i,j]*p[j]/den[i] d[j]<- sum( densum ) } } 下面是我算出來的 [1] 0.336093550 0.855872710 2.158927311 0.233847299 0.282627585 [6] 0.253739688 0.517757929 0.312419250 0.760612933 1.288101743 -- ※ 發信站: 批踢踢實業坊(ptt.cc), 來自: 219.91.75.186 (臺灣) ※ 文章網址: https://www.ptt.cc/bbs/R_Language/M.1574656447.A.C05.html ※ 編輯: empireisme (219.91.75.186 臺灣), 11/25/2019 12:39:29
→
11/25 13:44, 4年前 , 1F
sum(densum)不應該和for(i)同層
11/25 13:44, 1F
→
11/25 13:44, 4年前 , 2F
你先想一下j=1;i=1時densum是什麼就能明白有問題
11/25 13:44, 2F
→
11/25 14:09, 4年前 , 3F
11/25 14:09, 3F
→
11/25 14:11, 4年前 , 4F
其實它就只是每個i都取代一個d[j] 但是最後一個i
11/25 14:11, 4F
→
11/25 14:11, 4年前 , 5F
才會填一次正確值XD 是不影響結果拉...
11/25 14:11, 5F
→
11/25 14:12, 4年前 , 6F
這個問題用sweep / colSums / rowSums就能解決了
11/25 14:12, 6F
→
11/25 14:12, 4年前 , 7F
原po再自己看一下我上面貼的程式~ 有問題再發問
11/25 14:12, 7F
謝謝各位回答，我是沒用過sweep,只覺得這種迴圈寫的很煩，現在沒電腦，回家再試試 看
→
11/25 14:15, 4年前 , 8F
嗯後來發現了謝謝
11/25 14:15, 8F
所以我程式應該是沒寫錯的? 因為計算量太大也不好手算qq ※ 編輯: empireisme (101.12.41.101 臺灣), 11/25/2019 14:17:10 ※ 編輯: empireisme (101.12.41.101 臺灣), 11/25/2019 14:19:47
→
11/25 14:27, 4年前 , 9F
不算寫錯 但是就多了不少操作 是可以省掉的
11/25 14:27, 9F
→
11/25 14:28, 4年前 , 10F
這種情況建議拿 2 x 2的矩陣來驗算會好一點
11/25 14:28, 10F
→
11/25 14:29, 4年前 , 11F
sweep正如它字面上的意思就掃過去...
11/25 14:29, 11F
→
11/25 14:29, 4年前 , 12F
自己玩玩看吧
11/25 14:29, 12F
我算過小矩陣 但就是怕有萬一qq 比如說是不是對數學公式的理解錯誤 ※ 編輯: empireisme (101.12.41.101 臺灣), 11/25/2019 14:34:45
→
11/25 14:39, 4年前 , 13F
只看式子，確實多算很多。rowSums colSums先練練看
11/25 14:39, 13F
ok 看了c大的程式，看起來確實我沒有算錯但多算很多了 ※ 編輯: empireisme (101.12.41.101 臺灣), 11/26/2019 10:09:21
→
12/01 09:32, 4年前 , 14F
temp = t(A)*p
12/01 09:32, 14F
→
12/01 09:32, 4年前 , 15F
den = colSums(temp)
12/01 09:32, 15F
→
12/01 09:32, 4年前 , 16F
r2 = t(temp)/den
12/01 09:32, 16F
→
12/01 09:33, 4年前 , 17F
d = colSums(r2)
12/01 09:33, 17F
→
12/01 09:33, 4年前 , 18F
快一點，記憶體用量少一些
12/01 09:33, 18F
推
01/01 00:04, 5年前 , 19F
我的寫法是 simple=function(A,p){t=p/apply(A,1,function
01/01 00:04, 19F
→
01/01 00:04, 5年前 , 20F
(x)x%*%p);apply(A,2,function(x)x%*%t} 做 10000*10000矩
01/01 00:04, 20F
→
01/01 00:04, 5年前 , 21F
陣是4秒
01/01 00:04, 21F
‣ 更多 empireisme 的文章
文章代碼(AID): #1Tsrc_m5 (R_Language)
討論串 (同標題文章)
以下文章回應了本文：
完整討論串 (本文為第 1 之 2 篇)：
排序： 最舊先 | 最新先 | 留言數
在新視窗開啟完整討論串 (共2篇)
R_Language 近期熱門文章
更多 近期熱門文章 >>
PTT數位生活區 即時熱門文章
更多 即時熱門文章 >>
‣ 更多 empireisme 的文章
文章代碼(AID): #1Tsrc_m5 (R_Language)