qbeta function in R
It's here (search the sources for the appropriately named qbeta.c) but I can't guarantee you'll understand it easily: Martin Maechler (and others) works hard to make these the best in the business, but that performance comes at a price, here paid in opacity: http://svn.r-project.org/R/trunk/src/nmath/qbeta.c Michael
On Fri, Mar 9, 2012 at 3:08 PM, Anamika Chaudhuri <canamika at gmail.com> wrote:
Hi David: Thanks, I actually realized it was error on my part. I didnt think there was a problem with R-code, just wanted to get a better understanding of the function qbeta. Anamika On Fri, Mar 9, 2012 at 3:04 PM, David Winsemius <dwinsemius at comcast.net>wrote:
On Mar 9, 2012, at 2:48 PM, Anamika Chaudhuri wrote: ?HI All:
Does anyone know the code behind the qbeta function in R?
Well, yes, but don't you think it would be wise to question whether your code might be the problem rather than the R code? ?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])**)
Isn't this going to make x get longer with each pass through the loop? I think your parameter are then going to be interpreted as values of "x". Looks like "user error" to me. -- David
# 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]]
______________________________**________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/**listinfo/r-help<https://stat.ethz.ch/mailman/listinfo/r-help> PLEASE do read the posting guide http://www.R-project.org/** posting-guide.html <http://www.R-project.org/posting-guide.html> and provide commented, minimal, self-contained, reproducible code.
David Winsemius, MD West Hartford, CT
? ? ? ?[[alternative HTML version deleted]]
______________________________________________ R-help at 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.