Hello,
I am trying to use optim() on a function involving a summation. My
function basically is a thinned poisson likelihood. I have two parameters
and in most cases optim() does a fine job of getting the minima. I am
simulating my data based on pre specified parameters, so I know what I
should be getting. However when my true parameters fall in a particular
range, optim() gives incorrect results. I have generated a grid of
parameter values and calculated my likelihood through those to see that the
values that optim() gives is clearly not correct. I can see that my
likelihood does in fact have a unique maximum. Any ideas why this might be?
The data given below was generated such that the true parameters should be
(0.3001, -1.8971). Here is an example piece of data and the function:
#function to maximize
likN_alpha <- function(params,N){
thetas <- exp(params)
k <- length(thetas)
N <- c(rep(0,(k-1)),N)
l <- length(N)
lik <- 0
for(i in (k):(l-1)){
lambda <- thetas%*%N[i:(i-k+1)]
lik <- -lambda + N[i+1]*log(lambda) + lik
}
return(-lik)
}
# data to maximize over
N <- c( 3, 3, 10, 19, 36, 54, 78,116,177, 265, 388, 598,
890,1328,1910,2736,3982,5908,8471,12440,17964,26207,37688,54795,79270,114752,166594,242438,
352753,512054)
#optim() command
optim(log(1.5*rep(1/2,2)),likN.alpha,N=N)
# If I use constrained optimization and a slightly different
parameterization, then the results are fine, at least in this case, but not
always.
Thanks for any help you might be able to offer!
laura
______________________________________________
[email protected] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html