2017-01-09 305 views
0

積分を含む尤度を解くためにRでoptim()を使うのに問題があります。 L-BFGS-Bが有限の値fnを必要とするというエラーが表示されます(「par-c(0.1、0.1)、LLL、method = "L-BFGS-B" 「」以下は、私のコードは次のとおりです。。あなたの助けをRの最適化()の最適化は、 'fn'の有限値を必要とします。

s1=c(1384,1,1219,1597,2106,145,87,1535,290,1752,265,588,1188,160,745,237,479,39,99,56,1503,158,916,651,1064,166,635,19,553,51,79,155,85,1196,142,108,325 
    ,135,28,422,1032,1018,128,787,1704,307,854,6,896,902) 


LLL=function (par) { 

    integrand1 <- function(x){ (x-s1[i]+1)*dgamma(x, shape=par[1], rate=par[2]) } 
    integrand2 <- function(x){ (-x+s1[i]+1)*dgamma(x, shape=par[1],rate=par[2]) } 



    likelihood = vector() 

    for(i in 1:length(s1)) {likelihood[i] = 
    log(integrate(integrand1,lower=s1[i]-1,upper=s1[i])$value+ integrate(integrand2,lower=s1[i],upper=s1[i]+1)$value) 
    } 

    like= -sum(likelihood) 
    return(like) 

} 




optim(par=c(0.1,0.1),LLL,method="L-BFGS-B", lower=c(0,0)) 

おかげ

ベスト、

YM

+0

あなたの質問は読めないです。コードを表示する正しい方法を使用してください。 – Bhas

答えて

0

あなたが提供されたパラメータの下限で計算した目的関数が無限大であります。

LLL(c(0,0)) 
# [1] Inf 

だからこそ、L-BFGS-Bは失敗します。別の下限(たとえば、c(0.001,0.001))を試してみると、解決策が得られます。パラメータについて95%信頼区間を取得するには

optim(par=c(0.1,0.1),LLL,method="L-BFGS-B", lower=c(0.001,0.001)) 

$par 
[1] 0.6865841 0.0010000 

$value 
[1] 369.5532 

$counts 
function gradient 
     14  14 

$convergence 
[1] 0 

$message 
[1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH" 

はこれを試してみてください。

res <- optim(par=c(0.1,0.1),LLL,method="L-BFGS-B", lower=c(0.005,0.005), hessian=TRUE) 
n <- length(s1) 
res$par # solution 
# [1] 1.900928 0.005000 
res$par - 1.96*sqrt(diag(solve(res$hessian)))/n # lower limit for 95% confint 
# [1] 1.888152372 0.004963286 
res$par + 1.96*sqrt(diag(solve(res$hessian)))/n # upper limit for 95% confint 
# [1] 1.913703040 0.005036714 

はこの記事を参照してください。http://www.ms.uky.edu/~mai/sta321/MLEexample.pdf

+0

あなたは答えに私を打つ。私は全く同じソリューションを持っていました。 – Bhas

+0

@Bhas :) –

+0

@sandipanとBhasお返事いただきありがとうございます。 2ガンマパラメータの信頼区間を計算できるようにヘッセ行列をどのように出力できるかも分かりますか?私はoptim()で "hessian = TRUE"を追加しようとしました:optim(par = c(0.1,0.1)、LLL、method = "L-BFGS-B"、lower = c(0.001、: - 有限差分値[2] –