Skip to content
Prev 140070 / 398502 Next

problem with optim and integrate

kathie wrote:
Kathryn,

at first, please make your code readable and clarify it. I tried as follows:

fc <- function(theta, alpha){
   numint <- function(z, theta){
     (theta[1] * dnorm(z, mean = theta[2], sd = theta[3]) +
      theta[4] * dnorm(z, mean = theta[5], sd = theta[6]) +
      (1 - (theta[1] + theta[4])) *
           dnorm(z, mean = theta[7], sd = theta[8]))^alpha
   }
   integrate(numint, lower = -Inf, upper = Inf, theta = theta)$value
}

theta0 <- c(0.5, 0, sqrt(10), 0.25, -0.3, sqrt(0.05), 0.3, sqrt(0.05))

optim(theta0, fc, method = "L-BFGS-B", hessian = TRUE,
       lower = c(rep(c(0,-Inf,0), 2), -Inf, 0),
       upper = c(rep(c(1,Inf,4), 2), Inf, 4),
       alpha = 1.3)


Using some debugging tools (such as a simple browser() call or setting 
options(error=recover)), you will find that numint() can be NaN, because 
(1-(theta[1] + theta[4])) can be negative and hence the sum can be 
negative and a negative value to the power of alpha is not defined.

You may need to penalize for such cases.

Best,
Uwe Ligges