Гессен-матрица в optim() в R

У меня возникают проблемы с использованием метода optim() в R для решения задачи с вероятностью интеграла и получения матрицы hessian для двух параметров. Алгоритм сходится, но я получаю сообщение об ошибке, когда использую параметр hessian = TRUE в optim(). Ошибка:Гессен-матрица в optim() в R

Ошибка в интеграции (integrand1, низшие = S1 [I] - 1, верхний = S1 [I]): нефинитных значение функции

Также имел предупреждающее сообщение НАН

Вот мой код:

s1=c(1384,1,1219,1597,2106,145,87,1535,290,1752,265,588,1188,160,745,237,479,39,99,56,1503,158,916,651,1064,166,635,19,553,51,79,155,85,1196,142,108,325 
,135,28,422,1032,1018,128,787,1704,307,854,6,896,902) 

LLL=function (par) { 

    integrand1 <- function(x){ (x-s1[i]+1)*dgamma(x, shape=par[1], rate=par[2]) } 
    integrand2 <- function(x){ (-x+s1[i]+1)*dgamma(x, shape=par[1],rate=par[2]) } 



    likelihood = vector() 

    for(i in 1:length(s1)) {likelihood[i] = 
    log(integrate(integrand1,lower=s1[i]-1,upper=s1[i])$value+ integrate(integrand2,lower=s1[i],upper=s1[i]+1)$value) 
    } 

    like= -sum(likelihood) 
    return(like) 

} 




optim(par=c(0.1,0.1),LLL,method="L-BFGS-B", lower=c(0,0)) 
optim(par=c(0.1,0.1),LLL,method="L-BFGS-B", lower=c(0,0), hessian=TRUE)

Спасибо за вашу помощь!

источник

2017-01-09 Y. Ma

Я считаю, что вы оцениваете 'dgamma (0, shape = .1, scale = .1)', который является 'Inf', при интеграции для' i == 2' из 's1 [2] -1' (0) на 's1 [2]' (1). Что произойдет, если вы исключите '1' из' s1'? – Thales

@Thales Спасибо. –

Вам сказали вчера использовать 'lower = c (.001, .001)'. Почему вы не используете это для аргумента 'lower'? Измените 'lower' на' c (0.01,0.01) ', и вы увидите гессиан. – Bhas

optim минимизирует функцию. Вы можете построить функцию правдоподобия с аргументом rate. Ему нужно немного возиться, чтобы получить сюжет. Делают это так:

z2 <- function(rate) { 
    par <- numeric(2) 
    par[1] <- .68 
    par[2] <- rate 
    y <- LLL(par) 
    y 
} 

z1 <- Vectorize(z2,vectorize.args="rate") 

curve(z1, from=.001,to=1)

Вы увидите, что функция является минимальным для наименьшего значения для rate. То же самое, если вы меняете from на .1. Я не могу судить, действительна ли оценка.

источник

2017-01-09 19:53:53 Bhas

Еще раз спасибо за вашу большую помощь. –

Гессен-матрица в optim() в R

ответ

Смежные вопросы