基於牛頓拉夫森的最大似然估計和矩量法

Question

我正在進行一個小型模擬研究，以檢驗矩方法和漸近最大似然估計量的性質。

矩估計方法很容易獲得（它由我的代碼的第二行給出），但對於我必須編寫牛頓 - 拉夫遜算法（對於一個樣本來說工作得非常好）。因爲它以這種方式享有一些最佳統計特性，因此需要使用矩估計方法作爲起點（a0）。這裏去

x<-rbeta(500,0.5,3) 
mom<-3*mean(x)/(1-mean(x)) 

mlea<-function(x,a0,numstp=100,eps=0.001){ 
    n=length(x) 
    numfin=numstp 
    ic=0 
    istop=0 
    while(istop==0){ 
     ic=ic+1 
     lprime=n/a0+n/(a0+1)+n/(a0+2)+sum(log(x)) 
     ldprime=-n/a0^2-n/(a0+1)^2-n/(a0+2)^2 
     a1=a0-(lprime/ldprime) 
     check=abs((a1-a0)/a0) 
     if(check<eps){istop=1} 
     a0=a1 
    } 
    list(a1=a1,check=check,realnumstps=ic) 
} 

這適用於一個樣本，但我可以獲得這些估計說1000樣本？我怎麼能輕易地概括這個過程？我的主要困難是由於母親需要媽媽作爲起點，而這兩方面都需要從同一個樣本計算出來。

預先感謝您。

Answer 1

我想你想這樣做？

n<-100 
replicate(n, { 
    x<-rbeta(500,0.5,3) 
    mom<-3*mean(x)/(1-mean(x)) 
    mlea(x, mom) 
}) 

現在一個數字向量不會被返回，因爲你的函數mlea返回一個列表。比方說，從該列表中你真正關心的值是A1，那麼你可以做

n<-100 
replicate(n, { 
    x<-rbeta(500,0.5,3) 
    mom<-3*mean(x)/(1-mean(x)) 
    mlea(x, mom)$a1 
}) 

通知我稱之爲$ A1在函數調用結束。

所以，這是怎麼回事就在這裏被複制會拉從您的Beta分佈500次新的觀測中重複每次迭代（這將重複n次），然後計算依據是X的母親，然後給球菌mleA

的結果

我的結果？

replicate(3, { 
      x<-rbeta(500,0.5,3) 
      mom<-3*mean(x)/(1-mean(x)) 
      list(a1=mlea(x, mom)$a1, mom=mom) 
     }) 
# [,1]  [,2]  [,3]  
#a1 0.494497 0.522325 0.5153832 
#mom 0.4955767 0.5083678 0.5206997 

其中這裏的每一列是一個觀察