2011-07-04 51 views
-1

請問你能幫我糾正我的循環功能嗎?在R中優化循環:你能幫我調整我的循環代碼嗎?

我想優化估計al_j; au_j; sigma_j; b_j通過查看0到20,21到40,41到60個數據點。

最終結果應具有各估計中的4列和4行,每行的0〜20的,21〜40,41〜60

MY代碼是:

n=20 
runs=4 
out=matrix(0,nrow=runs) 

llik = function(x) 
{ 
    al_j=x[1]; au_j=x[2]; sigma_j=x[3]; b_j=x[4] 
    sum(na.rm=T, 
     ifelse(a$R_j< 0, -log(1/(2*pi*(sigma_j^2)))- 
         (1/(2*(sigma_j^2))*(a$R_j+al_j-b_j*a$R_m))^2, 
     ifelse(a$R_j>0 , -log(1/(2*pi*(sigma_j^2)))- 
         (1/(2*(sigma_j^2))*(a$R_j+au_j-b_j*a$R_m))^2, 
         -log(pnorm(au_j,mean=b_j*a$R_m,sd=sqrt(sigma_j^2))- 
         pnorm(au_j,mean=b_j*a$R_m,sd=sqrt(sigma_j^2))))) 

    ) 

} 

start.par = c(0, 0, 0.01, 1) 
out1 = optim(llik, par=start.par, method="Nelder-Mead") 


for (i in 1: runs) 
{ 
index_start=20*(i-1)+1 
index_end= 20*i 
out[i]=out1[index_start:index_end] 
} 
out 

我數據:

R_j   R_m 
-0.0625  0.002320654 
0    -0.004642807 
0.033333333 0.005936332 
0.032258065 0.001060848 
0    0.007114057 
0.015625  0.005581558 
0    0.002974794 
0.015384615 0.004215271 
0.060606061 0.005073116 
0.028571429 -0.006001279 
0   -0.002789594 
0.013888889 0.00770633 
0    0.000371663 
0.02739726 -0.004224228 
-0.04   0.008362539 
0   -0.010951605 
0    0.004682924 
0.013888889 0.011839993 
-0.01369863 0.004210383 
-0.027777778 -0.04658949 
0    0.00987272 
-0.057142857 -0.062203157 
-0.03030303 -0.119177639 
0.09375  0.077054642 
0   -0.022763619 
-0.057142857 0.050408775 
0    0.024706076 
-0.03030303 0.004043701 
0.0625  0.004951088 
0   -0.005968731 
0   -0.038292548 
0    0.013381097 
0.014705882 0.006424728 
-0.014492754 -0.020115626 
0    -0.004837891 
-0.029411765 -0.022054654 
0.03030303 0.008936428 
0.044117647 8.16925E-05 
0    -0.004827246 
-0.042253521 0.004653096 
-0.014705882 -0.004222151 
0.029850746 0.000107267 
-0.028985507 -0.001783206 
0.029850746 -0.006372981 
0.014492754 0.005492374 
-0.028571429 -0.009005846 
0    0.001031683 
0.044117647 0.002800551 
+0

這是作業嗎? –

+0

你能概括這個,以便它可能涉及更大的R社區? – Chase

+0

我已經創建了一個代碼,它使用R中的優化估計變量。我希望代碼執行一個循環(正如我在底部嘗試的那樣),但是我錯過了某些操作? – rder

回答

0

當我嘗試和運行代碼我得到錯誤信息「a未找到」當llik被調用。你應該重構你的代碼,這樣所有的變量llik或者被傳入它或者在函數中創建。

out[i]=out1[index_start:index_end]也看起來可疑,因爲您正在爲單個位置分配一個值的範圍。

+0

a是數據的地方。代碼如下:a = read.table(「E:/data.txt」,header = T)> attach(a)在a中是R_j和R_m。謝謝愛德華 – rder