[R] MLE with parameters restricted to a defined range using bbmle

Bernardo Santos bernardo_brandaum at yahoo.com.br
Wed Dec 10 14:35:15 CET 2014

```Thanks, Rolf and Ben,
Both solutions worked, but I finished by contraining parameters values using optim().
Best,Bernardo Niebuhr

Em Segunda-feira, 8 de Dezembro de 2014 18:36, Rolf Turner <r.turner at auckland.ac.nz> escreveu:

I know nothing about the bbmle package and its mle2() function, but it
is a general truth that if you need to constrain a parameter to be
positive in an optimisation procedure a simple and effective approach is
to reparameterize using exp().

I.e. represent xmin as exp(lxmin) (say) and use lxmin as the argument

This strategy rarely if ever fails to work.

cheers,

Rolf Turner

On 09/12/14 09:04, Bernardo Santos wrote:
> Dear all,
> I am fitting models to data with mle2 function of the bbmle package.In specific, I want to fit a power-law distribution model, as defined here (http://arxiv.org/pdf/cond-mat/0412004v3.pdf), to data.
> However, one of the parameters - xmin -, must be necessarily greater than zero. What can I do to restrict the possible values of a parameter that are passed to the optimizer?
> Here there is a sample of my code:
> library(bbmle)
>
> # Creating data
> set.seed(1234)
> data <- rexp(1000, rate = 0.1) # The fit will not be too good, but it is just to test
>
> # Creating the power-law distribution density function
> dpowlaw <- function(x, alfa, xmin, log=FALSE){
>    c <- (alfa-1)*xmin^(alfa-1)
>    if(log) ifelse(x < xmin, 0, log(c*x^(-alfa)))
>    else ifelse(x < xmin, 0, c*x^(-alfa))
> }
> # Testing the function
> integrate(dpowlaw, -Inf, Inf, alfa=2, xmin=1)
> curve(dpowlaw(x, alfa=2.5, xmin=10), from=0, to=100, log="")
> curve(dpowlaw(x, alfa=2.5, xmin=1), from=1, to=100, log="xy")
>
> # Negative log-likelihood function
> LLlevy <- function(mu, xmin){
>    -sum(dpowlaw(data, alfa=mu, xmin=xmin, log=T))
> }
>
> # Fitting model to data
> mlevy <- mle2(LLlevy, start=list(mu=2, xmin=1))
> The result of model fitting here is Coefficients:
>        mu      xmin
> -916.4043  890.4248
> but this does not make sense!xmin must be > 0, and mu must be > 1.What should I do?