Hi, R guys:
I'm using L-BFGS-B method of optim for minimization problem. My function
called besselI function which need non-negative parameter and the besselI
will overflow if the parameter is too large. So I set the constraint box
which is reasonable for my problem. But the point outside the box was
test, and I got error. My program and the error follows. This program
depends on CircStats package.
Anyone has any idea about this?
Thanks in advance.
Li
#################### source code ###################################
Dk2<- function(pars,theta)
{
kappa<- pars[1]; mu<- pars[2];
IoK<- besselI(kappa, nu=0);
res<- besselI(2*kappa, nu=0)/2/IoK^2 -
mean(exp(kappa*cos(theta-mu)))/IoK;
if(is.na(res)||is.infinite(res)){
print(pars);
# assign("Theta", theta, env=.GlobalEnv);
}
return(res);
}
mse.Dk2<- function(pars, s, n)
{
sum.est <- SSE <- numeric(2);
j<- 0;
while(j<=n){
theta<- rvm(s, pi, k=pars[1]) - pi;
est<- optim(par=pars, fn=Dk2, lower=c(0.001, -pi), upper=c(10,
pi), method="L-BFGS-B", theta=theta);
i<- 0;
while(est$convergence!=0 && i< 30){
est<- optim(par=est$par, fn=Dk2, lower=c(0.001, -pi),
upper=c(10, pi), method="L-BFGS-B", theta=theta);
i<- i+1;
}
if(est$convergence!=0) {
#print(j);
next;
}
else { j<- j+1; }
#est<- nlm(p=pars, f=Dk2, theta=theta);
mu.hat<- est$par[2];
while(mu.hat< -pi) mu.hat<- mu.hat + 2*pi;
while(mu.hat > pi) mu.hat<- mu.hat -2*pi;
est<- c(est$par[1], mu.hat);
sum.est <- sum.est + est;
SSE <- SSE + (est - pars)^2;
}
Est <- sum.est/n;
Bias<- Est - pars;
MSE<- SSE/n;
res<- c(Kappa=pars[1], Kappa.hat= Est[1], Kappa.Bias=Bias[1],
Kappa.MSE=MSE[1], Mu.hat=Est[2], Mu.MSE=MSE[2])
return(res);
}
kappas <- c(0.01, 0.05, 0.1, 0.20, 0.5, 1, 2, 5);
N<- 10000;
for ( s in c(5, 10, 20, 30, 50)){
cat("\nSample size = ", s);
cat("\n=====================================\n");
res<- NULL;
for(i in 1:8){
res<- rbind(res, mse.Dk2(c(kappas[i], 0), s, N));
}
print(round(res,4));
}
#Error message. -32.7 is far lower then the lower limit 0.001.
Sample size = 5
=====================================
[1] -32.736857 -3.141593
Error in optim(par = pars, fn = Dk2, lower = c(0.001, -pi), upper = c(10,
:
L-BFGS-B needs finite values of fn
In addition: Warning messages:
1: NaNs produced in: besselI(x, nu, 1 + as.logical(expon.scaled))
2: NaNs produced in: besselI(x, nu, 1 + as.logical(expon.scaled))
______________________________________________
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
