[R] Problem with fitdistr for beta

[Spencer Graves]

My standard work-around for the kind of problem you identified is to 
shrink the numbers just a little towards 0.5.  For example:

 > library(MASS)
 > a <- rbeta(100,0.1,0.1)
 > fitdistr(x=a, "beta", start=list(shape1=0.1,shape2=0.1))
Error in optim(start, mylogfn, x = x, hessian = TRUE, ...) :
	Function cannot be evaluated at initial parameters

 > r.a <- range(a)
 > c0 <- 0
 > c1 <- 1
 > if(r.a[1]==0)c0 <- min(a[a>0])
 > if(r.a[2]==1)c1 <- max(a[a<1])
 > c. <- c(c0, 1-c1)
 > if(any(c.>0))c. <- min(c.[c.>0])
 > c.
[1] 2.509104e-14
 > fitdistr(x=0.5*c.[1] + (1-c.[1])*a, "beta", 
start=list(shape1=0.1,shape2=0.1))
      shape1       shape2
   0.08728921   0.10403875
  (0.01044863) (0.01320188)
 >
hope this helps.  spencer graves

Prof Brian Ripley wrote:
> rbeta(100,0.1,0.1) is generating samples which contain 1, an impossible 
> value for a beta and hence the sample has an infinite log-likelihood.
> It is clearly documented on the help page that the range is 0 < x < 1.
> However, that is not so surprising as P(X > 1-1e-16) is about 1% and hence 
> values will get rounded to one.
> 
> The same would happen for a value of 0.
> 
> Your code is syntactically incorrect, at least as received here.
> 
> On Fri, 4 Jul 2003, Agustin Lobo wrote:
> 
> 
>>I have the following problem:
>>
>>I have a vector x of data (0<x<=1 ) with
>>a U-shaped histogram and try to fit a beta
>>distribution using  fitdistr. In fact,
>>hist(rbeta(100,0.1,0.1)) looks a lot like
>>my data.
>>
>>The equivalent to
>>the example in the manual
>>sometimes work:
>>
>>
>>>a <- rbeta(100,0.1,0.1)
>>>fitdistr(x=a, "beta", start=list(shape1=0.1,shape2=0.1))1)
>>>     shape1       shape2
>>
>>  0.09444627   0.12048753
>> (0.01120670) (0.01550129)
>>
>>but sometimes does not:
>>
>>>a <- rbeta(100,0.1,0.1)
>>>fitdistr(x=a, "beta", start=list(shape1=0.1,shape2=0.1))1)
>>>Error in optim(start, mylogfn, x = x, hessian = TRUE, ...) :
>>
>>        Function cannot be evaluated at initial parameters
>>
>>Unfortunately, my data fall in the second case
>>
>>I've searched for any weird value that be present in the
>>cases in which fitdistr exits with the error message, but
>>could not find any.
>>
>>Any help?
>>(please if anyone answers be sure to answer to my address as well,
>>I cannot subscribe to the list)
>>
>>Thanks
>>
>>Agus
>>
>>Dr. Agustin Lobo
>>Instituto de Ciencias de la Tierra (CSIC)
>>Lluis Sole Sabaris s/n
>>08028 Barcelona SPAIN
>>tel 34 93409 5410
>>fax 34 93411 0012
>>alobo at ija.csic.es
>>
>>______________________________________________
>>R-help at stat.math.ethz.ch mailing list
>>https://www.stat.math.ethz.ch/mailman/listinfo/r-help
>>
> 
>