[R-br] [Dúvida] Hessiana pela função optim

[Pedro Rafael]

Obrigado Wagner pela dica. Não conhecia o pacote bbmle. Na verdade sempre
fui tentado a sempre implementar tudo. Porém, seria uma boa tentar utilizar
o pacote.

Att.

Em 5 de agosto de 2016 09:39, Wagner Bonat <wbonat em gmail.com> escreveu:

> Esse problema pode ocorrer por muitos motivos incluindo apenas uma má
> aproximação do hessiano.
> Eu recomendo algumas coisas:
> 1 - Calcule a verossimilhança perfilhada de cada parâmetro. Isso vai de
> dar certeza que o algoritmo convergiu e te mostrar a cara da função. De
> acordo com o comportamento da profile likelihood vc pode tentar alguma
> reparametrização pra melhorar o comportamento do algoritmo numerico. Eu
> usaria o pacote bbmle que vai fazer isso automaticamente pra vc.
> 2 -  Calcule ao menos a primeira derivada da sua log-verossimilhança e use
> no BFGS.
> 3 - Calcule o hessiano numericamente porém separado, vc pode usar por
> exemplo o pacote numDeriv.
>
> Att
>
> Em 4 de agosto de 2016 22:39, Pedro Rafael via R-br <
> r-br em listas.c3sl.ufpr.br> escreveu:
>
>> Caros,
>>
>> sabemos que a matriz Hessiana multiplicada por -1 da função de
>> log-verossimilhança nos fornecem a matriz de informação observada que
>> converge assintoticamente para a matriz de informação esperada. A diagonal
>> principal da matriz inversa da informação esperada nos dá as variâncias os
>> estimadores de máxima verossimilhança.
>>
>> Um fato importante para o o questionamento que vou fazer é que estou
>> fazendo uso da matriz de informação observada de modo a ter ao menos uma
>> aproximação da variância dos estimadores de máxima verossimilhança obtidos
>> numericamente pelo método BFGS.
>>
>> O problema é que estou obtendo em alguns caros que a
>> solve(-diag(hessiana)) negativa o que não era de se esperar. Vejam se
>> possível o código.
>>
>> # PDF
>> pdf_ekww <- function(par,x){
>>   a = par[1]
>>   b = par[2]
>>   c = par[3]
>>   alpha = par[4]
>>   beta = par[5]
>>   g = dweibull(x = x, shape = alpha, scale = beta, log = FALSE)
>>   G = pweibull(q = x, shape = alpha, scale = beta, lower.tail = TRUE,
>> log.p = FALSE)
>>   a * b * c * g * G^(a-1) * (1-G^a)^(b-1) * (1-(1-G^a)^b)^(c-1)
>> }
>>
>> # Quantilica
>> sample_ekww <- function(n,par){
>>   a = par[1]
>>   b = par[2]
>>   c = par[3]
>>   alpha = par[4]
>>   beta = par[5]
>>   u = runif(n=n,min=0, max=1)
>>   qweibull(p = (1-(1-u^(1/c))^(1/b))^(1/a), shape = alpha, scale = beta,
>> lower.tail = TRUE, log.p = FALSE)
>> }
>> set.seed(1987)
>>
>> vector_par_true = c(1,1,1,1,1)
>> data = sample_ekww(n = 1000, par = vector_par_true)
>>
>> # Função de log-verossimilhança.
>> obj_ekww = function(par,x){
>>   sum(log(pdf_ekww(par,x)))
>> }
>>
>> # Maximizando log-verossimilhança.
>> result = optim(fn = obj_ekww, par = c(0.5,1.4,1,1,1), method = "BFGS", x
>> = data, hessian = TRUE, control=list("fnscale"=-1))
>>
>> diag(solve(-result$hessian))
>>
>> Observem que o primeiro elemento do vetor logo acima é negativo o que não
>> deveria ser verdade. Notem também que houve convergência segundo o critério
>> de parada da função optim, em que convergence é igual a zero. O que para
>> vocês podem estar provocando esse problema? A função objetivo é complicada
>> a ponto de provocar problemas na convergência do algoritmo vindo por sua
>> vez acarretar esse tipo de problema?
>>
>> Digo isso porque essas novas classes de distribuições de probabilidade
>> por apresentar diversos parâmetros produzem log-verossimilhanças muito
>> complicadas incluindo problemas de regiões aproximadamente planas bem como
>> problemas piores como log-verossimilhanças monótonas.
>>
>> Obrigado desde já,
>> Pedro Rafael.
>>
>>
>> _______________________________________________
>> R-br mailing list
>> R-br em listas.c3sl.ufpr.br
>> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
>> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça
>> código mínimo reproduzível.
>>
>
>
>
> --
> Wagner Hugo Bonat
> ------------------------------------------------------------
> ----------------------------------
> Department of Mathematics and Computer Science (IMADA)
> University of Southern Denmark (SDU) and
> Laboratório de Estatística e Geoinformação (LEG)
> Universidade Federal do Paraná (UFPR)
>
>
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20160808/1d5a59b2/attachment.html>