HI All:

Does anyone know the code behind the qbeta function in R?
I am using it to calculate exact confidence intervals and I am getting
'NaN' at places I shouldnt be. Heres the simple code I am using:

k<-3
> x<-NULL
> p<-rbeta(k,3,3)# so that the mean nausea rate is alpha/(alpha+beta)
> min<-10
> max<-60
> n<-as.integer(runif(3,min,max))
> for(i in 1:k)
+ x<-cbind(x,rbinom(5,n[i],p[i]))
>
> # Exact Confidence Interval
>
> l_cl_exact<-qbeta(.025, x, n-x+1)
Warning message:
In qbeta(p, shape1, shape2, lower.tail, log.p) : NaNs produced
> u_cl_exact<-qbeta(.975, x+1, n-x)
Warning message:
In qbeta(p, shape1, shape2, lower.tail, log.p) : NaNs produced
> x
     [,1] [,2] [,3]
[1,]    8   12   14
[2,]    5   15   13
[3,]    5   12   12
[4,]    8   21   12
[5,]    8   14   12
> n
[1] 10 36 31
> l_cl_exact
           [,1]      [,2]      [,3]
[1,] 0.44390454 0.2184996 0.2314244
[2,] 0.04667766       NaN 0.2454760
[3,] 0.05452433 0.1855618       NaN
[4,] 0.44390454 0.4862702 0.1855618
[5,] 0.10115053       NaN 0.2184996

Thanks for your help.
Anamika

	[[alternative HTML version deleted]]

