On 24/11/2007 7:43 AM, gallon li wrote:
> Suppose i want to compute a 95% highest density for a beta distribution
> beta(a,b)
>
> the two end points x1 and x2 shoudl satisfy the following two equations:
>
> pbeta(x1,a,b)-pbeta(x2,a,b)=95%
>
> dbeta(x1,a,b)=dbeta(x2,a,b)
>
> Is there any fast way to compute x1 and x2 in R?
I don't know if it's "fast", but uniroot can do it:
densitydiff <- function(lower, a, b, level=0.95) {
plower <- pbeta(lower, a, b)
pupper <- plower + level
upper <- qbeta(pupper, a, b)
return(dbeta(lower, a, b) - dbeta(upper, a, b))
}
a <- 2
b <- 10
x2 <- uniroot(densitydiff, c(0, qbeta(0.05, a,b)), a=a, b=b)$root
(and then calculate the upper limit the same way densitydiff did).
Duncan Murdoch
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
