[R] fitdistr, mle's and gamma distribution
Lourens Olivier Walters
lwalters at cs.uct.ac.za
Thu Oct 9 10:31:55 CEST 2003
Thanks, again I managed to get convergence, although looking at the
resultant Pareto distribution plotted in R, it doesn't seem that the log
likelihood estimated parameters are more accurate than the methods of
moments estimator parameters - I will have to look into it. For the
record, here is the code I used:
###########
data.sorted <- sort(data)
alpha.val <- data.sorted[1]
beta.val <- 1/((1/n) * sum(log(data.sorted/alpha.val)))
log.alpha.val <- log(alpha.val)
# Optimize log.del = log(0.0000001), which is the log(difference
# between alpha and log likelihood estimated alpha), chose small
# difference to start of with
log.del = -16.11809565
paretoLoglik <- function(params, negative=TRUE){
alpha <- data.sorted[1]+exp(params[1])
lglk <- sum(pareto.density(data.sorted, alpha=exp(alpha),
+ beta=exp(params[2]), log=TRUE))
if(negative)
return(-lglk)
else
return(lglk)
}
optim.list <- optim(c(log.del, log.beta.val), paretoLoglik,
+ method="BFGS")
pareto.param1 <- exp(optim.list$par[1]) + alpha.val
pareto.param2 <- exp(optim.list$par[2])
##########
On Wed, 2003-10-08 at 15:07, Spencer Graves wrote:
> I'm sorry, but I don't have time to read all your code.
However,
> I saw that you tested for x > alpha in your Pareto distribution
> example. Have you considered reparameterizing to estimate log.del =
> log(alpha-min(x))? Pass log.del as part of the vector of parameters
to
> estimate, then immediately inside the function, compute
>
> alpha <- (min(x)+exp(log.del))
>
> I've fixed many convergence problems with simple tricks like
this.
>
> hope this helps. spencer graves
More information about the R-help
mailing list