[R] Problem with fitdistr for beta
Spencer Graves
spencer.graves at pdf.com
Fri Jul 4 19:09:09 CEST 2003
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
>>
>
>
More information about the R-help
mailing list