Pedro,
A princípio, acho que você poderia tentar reparametrizar o modelo, colocando os
parâmetros em escala logarítmica (p. ex. a = exp(par[1]); b = exp(par[2]);
...), e tentar alguma outra função de otimização
(http://cran.r-project.org/web/views/Optimization.html).
Além disso, seria uma boa tentar simular dados desta distribuição com
parâmetros conhecidos e ver como se comportam as estimativas, em vez de tentar
ajustar um conjunto de dados diretamente.
Att.,
Rubem
________________________________
De: Pedro Rafael <[email protected]>
Para: [email protected]
Enviadas: Quarta-feira, 5 de Junho de 2013 14:30
Assunto: [R-br] [Dúvida] Estimação por Máxima Verossimilhança.
Senhores tenho uma dúvida. Na verdade não é dúvida, apenas quero sugestões.
Tenho uma função densidade de probabilidade "complicada". Trate-se de uma
distribuição chamada Kwmarashwamy Weibull Poisson. Tenho alguns bancos de dados
e gostaria de verificar o ajustamento dessa distribuição à estes bancos de
dados. Estou estimando os parâmetros pelo método de máxima verossimilhança.
Essa distribuição tem suporte nos reais positivos (x>0) e todos os seus
parâmetros são positivos. Optei em utilizar o método L-BFGS-G para restringir a
busca nos reais positivos. Segue abaixo o comando. Nesse exemplo não houve
convergência. Percebi que os chutes iniciais influenciam muito as estimativas
dos parâmetros no caso em que há convergência. Usando métodos de maximização
diferentes em muitos casos há grandes diferenças nas estimativas. Existe
alguma forma mais tranquila e direta para encontrar as estimativas pelo método
de máxima
verossimilhança em R? O código segue abaixo:
vero <- function(par,x){
a = par[1]
b = par[2]
c = par[3]
lambda = par[4]
beta = par[5]
-sum(log((a*b*c*lambda*(beta^c)*(x^(c-1))*((1-exp(-(x*beta)^c))^(a-1)) *
((1-(1-exp(-(beta*x)^c))^a)^(b-1)) *
exp(-lambda*(1-(1-(1-exp(-(beta*x)^c))^a)^b) -
(beta*x)^c))/(1-exp(-lambda))))
}
dados = c(17.23, 28.92, 33.00, 41.52,
42.12, 45.60, 48.80, 51.84, 51.96, 54.12, 55.56, 67.80, 68.64,
68.64,68.88,
84.12, 93.12, 98.64, 105.12, 105.84, 127.92, 128.04, 173.40)
optim(par=c(1,1,1,1,1),fn=vero,
method="L-BFGS-B",x=dados/1000,
lower=c(0.001,0.001,0.001,0.001,0.001), upper=c(Inf,Inf,Inf,Inf,Inf))
O engraçado é que essa distribuição é encaixada com a distribuição Weibull.
Quando tento ajustar a Weibull à esses dados há convergência. O que vocês acham
que devo fazer?
[ ],
Pedro Rafael Diniz Marinho.
_______________________________________________
R-br mailing list
[email protected]
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._______________________________________________
R-br mailing list
[email protected]
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.