Não consigo estimar por máxima verossimilhança (problemas com logaritmo de argumento negativo)

Caros colegas, estou com um problema que não consigo resolver... Na verdade acredito que já identifiquei o problemas, mas não sei como contornar esse tipo de situação no R. É o seguinte estou maximizando a função de log-verossimilhança de uma dist. de probabilidade, com o intuito de estimar os parâmetros da dist. Entretanto a otimização não está convergindo para o parâmetro conhecido. Essa dist. possui dois parâmetros, onde um deles o "z" assume valores menores que 2, e o outro parâmetro "eta" assume valores maiores que 0. Consegui simular situações para os casos que tenho 1<z<2 e para 0<z<1, entretanto para valores negativos de z, não estou conseguindo convergência para o parâmetro. Estou usando a função optim() do R, selecionando nessa função o método de Nelder-Mead (já usei outros métodos e nenhum convergiu). O problema que ocorre é que para z negativo, durante as iterações ocorrem o cálculo de logaritmos cujos argumentos estão com valores negativos, tenho quase certeza que o problema da não convergência está aí, aí eu gostaria de saber se tem como colocar alguma restrição durante o processo de iteração para que os argumentos das funções logarítmicas nunca sejam valores negativos. A seguir coloco o código que estou usando: ############################################## ##### Fornecendo os parâmetros ############################## n<-50 z<--1.95 eta<-0.1 z.est=vector() eta.est=vector() p=1 while (p<1001) { ###### Gerando valores da varíavel cujos parâmetros são z e eta u<-runif(n) y = eta*((1-(u)^((1-z)/(2-z)))/(1-z)) #### Obtendo as estimativas para z e eta ######################################## vero <- function(par,x){ z = par[1] eta = par[2] saida<-((sum(log(1 - ((1 - z)*x*(1/eta)))))/(1 - z)) + ((n)*log(2 - z)) + (n*log(1/eta)) return(-saida) } saida1<-optim(par=c(-1,4),fn=vero, method="Nelder-Mead",x=y ) z.est[p]<-saida1[1]$par[1] eta.est[p]<-saida1[1]$par[2] p=p+1 } mean(z.est) mean(eta.est) ############################################## Como estou replicando o procedimento mil vezes, irei utilizar a média das 1000 estimativas encontradas para z e para eta... Então fazendo isso tenho a seguinte saída:
mean(z.est) [1] -30177882
mean(eta.est) [1] 964132.7
Conforme observamos as estimativas estão muito distantes dos parâmetros reais definidos no início do código (z = -1.95 e eta = 0.1). Dessa forma, desconfio que o problema seja decorrente dos vários erros de logaritmos, que se repetem durante a interação e são da seguinte forma: 1: In log(1 - ((1 - z) * x * (1/eta))) : NaNs produzidos Então tem alguma forma de controlar esse problema que ocorre durante a interação??? Tenho que fazer de alguma forma com que os argumentos fiquem sempre positivos, ou seja: (1-((1-z)*x*(1/eta))) tem que ser sempre positivo... Por favor, caso vcs possam me ajudar, eu agradeço imensamente. Att, Romero.

Em problemas de otimização o primeiro conselho é ter uma função objetivo com argumentos no domínio dos reais. Tente colocar reparametrizações na sua função objetivo para garantir isso. Por exemplo, no caso da distribuição normal a função de log-verrossimilhança tem que ser escrita com mu=mu e sigma=exp(theta), não se estima diretamente sigma pois sigma>0 mas exp(theta)>0 para todos theta nos reais. A mesma idéia serve para binomial, onde p está no (0,1), não se otimiza em p mas em theta sendo que p=1/(1+exp(theta)) que vem da função de ligação log(p/(1-p))=theta. À disposição. Walmes.

Romero, faz tempo que não tenho esse tipo de problema na minha vida. A sugestão do Walmes me parece a mais correta. Outra abordagem que me parece ser razoável e tentar maximizar sem a função de log, pois deveria obter o mesmo resultado. Ou também usar uma abordagem bayesiana. Determinar uma distribuição a posteriori e então determinar as distribuições condicionais completas para os seus parâmetros. Sem dúvida esse tipo de approach vai te dar muito mais trabalho. Abs Vinicius Em 17 de novembro de 2014 20:47, walmes . <walmeszeviani@gmail.com> escreveu:
Em problemas de otimização o primeiro conselho é ter uma função objetivo com argumentos no domínio dos reais. Tente colocar reparametrizações na sua função objetivo para garantir isso. Por exemplo, no caso da distribuição normal a função de log-verrossimilhança tem que ser escrita com mu=mu e sigma=exp(theta), não se estima diretamente sigma pois sigma>0 mas exp(theta)>0 para todos theta nos reais. A mesma idéia serve para binomial, onde p está no (0,1), não se otimiza em p mas em theta sendo que p=1/(1+exp(theta)) que vem da função de ligação log(p/(1-p))=theta.
À disposição. Walmes.
_______________________________________________ R-br mailing list R-br@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.
-- *Vinicius Brito Rocha.* *Estatístico e Atuário * *M.Sc. Engenharia de Produção/PO* *Ph.D. Estudante Sistemas Computacionais*

Prezados, mais uma vez peço desculpas... Não consigo responder diretamente ao tópico dando o reply na mensagem que chega até a minha caixa de e-mail. Será que estou fazendo algum procedimento errado??? Como preciso muito solucionar meu caso, estou reenviando uma nova mensagem com o mesmo título do tópico anteriormente criado. Se existir a possibilidade de mesclar esse tópico com os outros dois criados eu ficaria agradecido. Se alguém puder me orientar tbm em como responder uma mensagem no R-BR eu agradeço, pois já participei de alguns tópicos aqui e acredito que sempre respondia ao e-mail que chegava na minha conta... Neste tópico vou colocar as mensagens postadas anteriormente para pelo menos seguir uma linha de raciocínio, e novamente peço para me orientarem como responder nesse tópico atual. Obrigado. ########################################### A seguir peço que leiam a situação atual do problema: ########################################### Caros colegas, estou com um problema que não consigo resolver... Na verdade acredito que já identifiquei o problemas, mas não sei como contornar esse tipo de situação no R. É o seguinte estou maximizando a função de log-verossimilhança de uma dist. de probabilidade, com o intuito de estimar os parâmetros da dist. Entretanto a otimização não está convergindo para o parâmetro conhecido. Essa dist. possui dois parâmetros, onde um deles o "z" assume valores menores que 2, e o outro parâmetro "eta" assume valores maiores que 0. Consegui simular situações para os casos que tenho 1<z<2 e para 0<z<1, entretanto para valores negativos de z, não estou conseguindo convergência para o parâmetro. Estou usando a função optim() do R, selecionando nessa função o método de Nelder-Mead (já usei outros métodos e nenhum convergiu). O problema que ocorre é que para z negativo, durante as iterações ocorrem o cálculo de logaritmos cujos argumentos estão com valores negativos, tenho quase certeza que o problema da não convergência está aí, aí eu gostaria de saber se tem como colocar alguma restrição durante o processo de iteração para que os argumentos das funções logarítmicas nunca sejam valores negativos. A seguir coloco o código que estou usando: ############################################## ##### Fornecendo os parâmetros ############################## n<-50 z<--1.95 eta<-0.1 z.est=vector() eta.est=vector() p=1 while (p<1001) { ###### Gerando valores da varíavel cujos parâmetros são z e eta u<-runif(n) y = eta*((1-(u)^((1-z)/(2-z)))/(1-z)) #### Obtendo as estimativas para z e eta ######################################## vero <- function(par,x){ z = par[1] eta = par[2] saida<-((sum(log(1 - ((1 - z)*x*(1/eta)))))/(1 - z)) + ((n)*log(2 - z)) + (n*log(1/eta)) return(-saida) } saida1<-optim(par=c(-1,4),fn=vero, method="Nelder-Mead",x=y ) z.est[p]<-saida1[1]$par[1] eta.est[p]<-saida1[1]$par[2] p=p+1 } mean(z.est) mean(eta.est) ############################################## Como estou replicando o procedimento mil vezes, irei utilizar a média das 1000 estimativas encontradas para z e para eta... Então fazendo isso tenho a seguinte saída:
mean(z.est) [1] -30177882
mean(eta.est) [1] 964132.7
Conforme observamos as estimativas estão muito distantes dos parâmetros reais definidos no início do código (z = -1.95 e eta = 0.1). Dessa forma, desconfio que o problema seja decorrente dos vários erros de logaritmos, que se repetem durante a interação e são da seguinte forma: 1: In log(1 - ((1 - z) * x * (1/eta))) : NaNs produzidos Então tem alguma forma de controlar esse problema que ocorre durante a interação??? Tenho que fazer de alguma forma com que os argumentos fiquem sempre positivos, ou seja: (1-((1-z)*x*(1/eta))) tem que ser sempre positivo... Por favor, caso vcs possam me ajudar, eu agradeço imensamente. Att, Romero. ########################################### ########################################### Mensagem anterior enviada por Wagner bonat: ************************************** Uma outra opcao é vc restringir o espaco de busca da optim() usando o metodo L-BFGS-B. Vc precisa identificar para quais valores do seu parametro o erro acontece e se certificar que o algoritmo nao caminhe nesta direcao. Uma segunda opcao é se for possivel calcule pelo menos a primeira derivada da funcao de log-verossimilhanca e a forneca ao algoritmo BFGS ou L-BFGS-B isso ajuda o algoritmo ficar sempre no espaco parametrico, porque a derivada nao deixa ele tentar valores meio que absurdos. Para uma implementacao definitiva eu aconselharia usar um Fisher-score muito eficiente e normalmente partindo de valores iniciais razoaveis vc sempre converge. 'E claro que para isso vc precisa calcular a matriz de derivadas segundas e tomar a esperanca, dependendo da situacao pode ser razoavel. ########################################### ########################################### Mensagem anterior enviada por Romero (rlmsf) ************************************** Olá Walmes e Vinícius, obrigado pelas dicas, pelo menos já estou buscando encontrar alguma forma de resolver, pois confesso, que já não estava vendo mais como sair do dilema. Vinicius, já tinha tentado fazer a otimização da função de verossimilhança (sem o log), entretanto, não sei o que acontece, mas as estimativas ficam todas iguais aos pontos iniciais, ou seja, para as mil replicações que estou fazendo, todos os resultados são iguais aos pontos iniciais... Com certeza, fazendo dessa forma, não deve estar chegando nem na segunda iteração (não sei o motivo). Walmes, a sua dica é boa, alguém já tinha falado em algo do tipo, mas eu não tinha entendido bem... acho que agora com sua explicação ficou mais claro, resumindo, por favor veja se é isso: preciso encontrar funções cujos possíveis resultado sejam pertencentes ao domínio dos meus dois parâmetros de interesse, com a exigência de que o argumento dessas funções sejam sempre valores pertencentes aos reais, certo?? Pensando assim, como o parâmetro z é sempre menor que 2 e o parâmetro eta é sempre um número positivo, então teoricamente posso usar a seguinte transformação: z= 2-exp(x) (argumento x real e resultado da função sempre menor que 2) e eta=exp (k) (argumento k real e resultado da função sempre maior que zero). É isso mesmo, vcs concordam??? Obrigado pela atenção, Romero. ########################################### ########################################### Mensagem anterior enviada por Vinicius Brito ************************************** Romero, faz tempo que não tenho esse tipo de problema na minha vida. A sugestão do Walmes me parece a mais correta. Outra abordagem que me parece ser razoável e tentar maximizar sem a função de log, pois deveria obter o mesmo resultado. Ou também usar uma abordagem bayesiana. Determinar uma distribuição a posteriori e então determinar as distribuições condicionais completas para os seus parâmetros. Sem dúvida esse tipo de approach vai te dar muito mais trabalho. Abs Vinicius ########################################### ########################################### Mensagem anterior enviada por Walmes ************************************** Em problemas de otimização o primeiro conselho é ter uma função objetivo com argumentos no domínio dos reais. Tente colocar reparametrizações na sua função objetivo para garantir isso. Por exemplo, no caso da distribuição normal a função de log-verrossimilhança tem que ser escrita com mu=mu e sigma=exp(theta), não se estima diretamente sigma pois sigma>0 mas exp(theta)>0 para todos theta nos reais. A mesma idéia serve para binomial, onde p está no (0,1), não se otimiza em p mas em theta sendo que p=1/(1+exp(theta)) que vem da função de ligação log(p/(1-p))=theta. À disposição. Walmes. ################################################### ################################################### Mensagem anterior enviada por Romero (rlmsf) (1º Tópico) ********************************************* Caros colegas, estou com um problema que não consigo resolver... Na verdade acredito que já identifiquei o problemas, mas não sei como contornar esse tipo de situação no R. É o seguinte estou maximizando a função de log-verossimilhança de uma dist. de probabilidade, com o intuito de estimar os parâmetros da dist. Entretanto a otimização não está convergindo para o parâmetro conhecido. Essa dist. possui dois parâmetros, onde um deles o "z" assume valores menores que 2, e o outro parâmetro "eta" assume valores maiores que 0. Consegui simular situações para os casos que tenho 1<z<2 e para 0<z<1, entretanto para valores negativos de z, não estou conseguindo convergência para o parâmetro. Estou usando a função optim() do R, selecionando nessa função o método de Nelder-Mead (já usei outros métodos e nenhum convergiu). O problema que ocorre é que para z negativo, durante as iterações ocorrem o cálculo de logaritmos cujos argumentos estão com valores negativos, tenho quase certeza que o problema da não convergência está aí, aí eu gostaria de saber se tem como colocar alguma restrição durante o processo de iteração para que os argumentos das funções logarítmicas nunca sejam valores negativos. A seguir coloco o código que estou usando: ############################################## ##### Fornecendo os parâmetros ############################## n<-50 z<--1.95 eta<-0.1 z.est=vector() eta.est=vector() p=1 while (p<1001) { ###### Gerando valores da varíavel cujos parâmetros são z e eta u<-runif(n) y = eta*((1-(u)^((1-z)/(2-z)))/(1-z)) #### Obtendo as estimativas para z e eta ######################################## vero <- function(par,x){ z = par[1] eta = par[2] saida<-((sum(log(1 - ((1 - z)*x*(1/eta)))))/(1 - z)) + ((n)*log(2 - z)) + (n*log(1/eta)) return(-saida) } saida1<-optim(par=c(-1,4),fn=vero, method="Nelder-Mead",x=y ) z.est[p]<-saida1[1]$par[1] eta.est[p]<-saida1[1]$par[2] p=p+1 } mean(z.est) mean(eta.est) ############################################## Como estou replicando o procedimento mil vezes, irei utilizar a média das 1000 estimativas encontradas para z e para eta... Então fazendo isso tenho a seguinte saída:
mean(z.est) [1] -30177882
mean(eta.est) [1] 964132.7
Conforme observamos as estimativas estão muito distantes dos parâmetros reais definidos no início do código (z = -1.95 e eta = 0.1). Dessa forma, desconfio que o problema seja decorrente dos vários erros de logaritmos, que se repetem durante a interação e são da seguinte forma: 1: In log(1 - ((1 - z) * x * (1/eta))) : NaNs produzidos Então tem alguma forma de controlar esse problema que ocorre durante a interação??? Tenho que fazer de alguma forma com que os argumentos fiquem sempre positivos, ou seja: (1-((1-z)*x*(1/eta))) tem que ser sempre positivo... Por favor, caso vcs possam me ajudar, eu agradeço imensamente. Att, Romero.
participantes (3)
-
Romero Luiz M. Sales Filho
-
Vinicius Brito Rocha
-
walmes .