Оценка среднего и стандартного отклонения нормального распределения в R.problem с функцией «optim»

Question

Я пытаюсь написать код в r, чтобы найти значения максимального правдоподобия (не лог-правдоподобия) для одномерное нормальное распределение. Я знаю, что существуют другие методы, но для дальнейших работ мне требуется глубокое понимание численной оптимизации. Когда я вызываю функцию «optim», кажется, что она не выполняет итерации вообще и возвращает значения, которые я передал в качестве начальных параметров. Этого не происходит, если я передаю оптимизатору функцию, которая вместо этого вычисляет лог-правдоподобие. Любая идея почему? Я не вижу, где моя ошибка. Я мог бы просто сказать, что продукт плотностей слишком близок к нулю, и калькулятор не может справиться с этим. Вот мой код. Спасибо большое!

set.seed(123); 
x=rnorm(10, mean = 2, sd = 5) 

LikeNormUnivar<-function(param,data){ 
    mu=param[1]; 
    sdev=param[2]; 
    densityvector=dnorm(data, mean = mu, sd = sdev, log = FALSE) 
    like=prod(densityvector) 
    return(-like) 
} 

theta.start = c(2,4) 
ans = optim(par=theta.start, fn=LikeNormUnivar, data=x,control=list(trace=TRUE), 
      method="BFGS") 
ans$par 

Answer 1

я понял, что происходит, добавив строку cat(mu,sdev,like,"\n") в соответствующем месте в функции, чтобы увидеть, что происходит. В основном, по шкале, на которой BFGS оценивает производные по конечным разностям, изменений недостаточно. Установка retol=1e-16 в списке control работает. Лучше попробуйте свести к минимуму отрицательный журнал правдоподобие ... и установка на шкале стандартного отклонения журнала также является хорошей идеей, например.

LikeNormUnivar <- function(param,data){ 
    mu=param[1] 
    sdev=exp(param[2]) 
    loglik=dnorm(data, mean = mu, sd = sdev, log = TRUE) 
    return(-sum(loglik)) 
}