问题描述
我有一个想要最大化的简单似然函数(来自均值= 0的正态分布).乐观主义者不断给我这个错误:optim中的错误(par = phi,fn = loglike,估计=估计,NULL,hessian = TRUE ,:非有限有限差分值[1]
I have a simple likelihood function (from a normal dist with mean=0) that I want to maximize. optim keeps giving me this error:Error in optim(par = p fn = loglike, estimates = estimates, NULL, hessian = TRUE, : non-finite finite-difference value [1]
这是我的数据和似然函数:
Here is my data and likelihood function:
y = [ -0.01472 0.03942 0.03592 0.02776 -0.00090 ]
C = a varcov matrix:
1.66e-03 -0.000120 -6.780e-06 0.000102 -4.000e-05
-1.20e-04 0.001387 7.900e-05 -0.000140 -8.000e-05
-6.78e-06 0.000079 1.416e-03 -0.000070 8.761e-06
1.02e-04 -0.000140 -7.000e-05 0.001339 -6.000e-05
-4.00e-05 -0.000080 8.761e-06 -0.000060 1.291e-03
我的对数似然函数是:lglkl = -.5 *(log(det(v))+(t(y)%%vi%%y)))其中v = phi * I + C和vi = inverse(v )和I = 5 * 5身份矩阵.
my log likelihood function is:lglkl = -.5*(log(det(v)) + (t(y)%%vi%%y))` where v = phi*I + C and vi=inverse(v) and I= 5*5 Identity matrix.
我正在尝试获取"phi"的最大似然估计.我以为这将是一个简单的优化问题,但仍然很麻烦.非常感谢您的帮助.提前致谢.我的代码如下:
I am trying to get the mle estimate for "phi". I thought this would be a simple optimization problem but am struggling. Would really appreciate any help. Thanks in advance. My code is below:
loglike <- function(phi,y) {
v = phi*I + C
vi = solve(v)
loglike = -.5*(log(det(v)) + (t(y)%*%vi%*%y))
return(-loglike)
}
phi = 0
parm <- optim(par=phi,fn=loglike,y=y,NULL,hessian = TRUE, method="L-BFGS-B",lower=0,upper=1000)
推荐答案
您遇到的错误是因为ϕ在经过一定数量的迭代后变为负(这表明算法未正确应用约束).同样,解决方案不会收敛到单个值,而是会在达到更新的协方差矩阵不再是正定值的情况之前跳到几个小值之间.在那个阶段,您得到det(v)< 0,并且log [det(v)]未定义. optim
算法在该阶段失败了.
The error you ran into is because ϕ becomes negative beyond a certain number of iterations (which indicates that the constraints are not being applied correctly by the algorithm). Also, the solution does not converge to a single value but jumps between a few small values before reaching a situation where the updated covariance matrix is no-longer positive definite. At that stage you get det(v) < 0 and log[det(v)] is undefined. The optim
algorithm bails out at that stage.
要查看发生了什么,请在下面的代码中使用maxit
和ndeps
参数.
To see what's happening, play with the maxit
and ndeps
parameters in the code below.
require("matrixcalc")
#-------------------------------------------------
# Log-likelihood function
#-------------------------------------------------
loglike <- function(phi, y) {
# Shift the covariance matrix
print(paste("phi = ", phi))
#v = phi*I + (1 - phi)*C
v = phi*I + C
stopifnot(is.positive.definite(v))
# Invert shifted matrix
vi = solve(v)
# Compute log likelihood
loglike = -.5*(log(det(v)) + (t(y) %*% vi %*% y))
print(paste("L = ", loglike))
return(-loglike)
}
#-------------------------------------------------
# Data
#-------------------------------------------------
y = c(-0.01472, 0.03942, 0.03592, 0.02776, -9e-04)
C = structure(c(0.00166, -0.00012, -6.78e-06, 0.000102, -4e-05, -0.00012,
0.001387, 7.9e-05, -0.00014, -8e-05, -6.78e-06, 7.9e-05,
0.001416, -7e-05, 8.761e-06, 0.000102, -0.00014, -7e-05,
0.001339, -6e-05, -4e-05, -8e-05, 8.761e-06, -6e-05, 0.001291),
.Dim = c(5L, 5L ))
#--------
# Initial parameter
#--------
I = diag(5)
phi = 50
#--------
# Minimize
#--------
parm <- optim(par = phi, fn = loglike, y = y, NULL, hessian = TRUE,
method = "L-BFGS-B", lower = 0.0001, upper = 1000,
control = list(trace = 3,
maxit = 1000,
ndeps = 1e-4) )
这篇关于r中的优化:非有限有限差分误差的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持!