Eu preciso parametrizar uma distribuição Weibull para alguns dados. Portanto, eu uso a Estimativa de Máxima Verossimilhança (MLE) do pacote fitdistrplus em R. No entanto, eu queria entender o que é feito no pacote, portanto, além de usar o pacote, tentei duas soluções manuais para verificar as propriedades do MLE fornecidas por fitdist.
Resumindo, minhas abordagens são:
(i) Use a função fitdist com o método "MLE"
(ii) Resolva as derivadas parciais da função de probabilidade
(iii) Minimize a probabilidade negativa usando a função otim
Primeiro, simule alguns dados:
n <- 1e4
set.seed(1)
dat <- rweibull(n, shape=0.8, scale=1.2)
Abordagem 1: aplique o pacote fitdistrplus:
library(fitdistrplus)
A1 <- fitdist(dat, "weibull", method="mle")$estimate
A1
shape scale
0.7914886 1.2032989
Abordagem 2:
Tendo como densidade Weibull
,
os derivados parciais são:
Pesquise as raízes das derivadas parciais acima:
weib1 <- function(c) { 1/c - sum(dat^c*log(dat))/sum(dat^c) + 1/n*sum(log(dat)) }
shape <- uniroot(weib1, c(0,10), tol=1e-12)$root
scale <- (1/n*sum(dat^shape))^(1/shape)
A2 <- c(shape, scale)
A2
[1] 0.7914318 1.2033179
Abordagem 3: procure os parâmetros que minimizam a probabilidade de log negativa:
fobj <- function(params){
-sum(log(dweibull(dat, params[1], params[2])))
}
A3 <- optim(c(0.5, 1), fobj)$par
A3
[1] 0.7913756 1.2032748
Comparando as abordagens, as estimativas dos parâmetros (A1, A2, A3) diferem na quarta casa decimal. Considerando a documentação do fitdist, eu esperava que A1 e A3 produzissem as mesmas estimativas, pois ambos usam otim.
Portanto, minhas perguntas são:
Qual é a função objetivo usada pelo fitdist e como eu poderia mudar a abordagem 3 para produzir exatamente as mesmas estimativas que o fitdist? E, em geral, qual seria a abordagem preferida, suponho que resolver as derivadas parciais é a abordagem mais limpa?
fonte