It's a pretty extreme case, try e.g. curve(pbeta(x, shape1, shape2), n=10001), and (probably -- I can't be bothered to work out the relation between beta shapes and F df parameters this morning...) outside what is normally encountered in statistical analyses. Notice though, that you have lower=FALSE in your Qbeta0, so you should have it in qbeta as well:
> qbeta(t,shape1,shape2, lower=FALSE) [1] 0.9999949 Warning message: In qbeta(t, shape1, shape2, lower = FALSE) : full precision may not have been achieved in 'qbeta' which of course is still wrong (I dont't think there is a problem in the other tail, Qbeta0(t,...., lower=TRUE) returns zero. You can see the effect also with curve(qbeta(x, shape1, shape2), n=10001, from=.99, to=1) which kind of suggests one of those regime-switching bugs, where different methods are used for various parts of the domain, and the cross-over is not done quite right. At any rate, qbeta is one of R's very basic workhorses, so we do want it to work right, so by all means file a report. -pd > On 26 Mar 2020, at 02:09 , Ben Bolker <bbol...@gmail.com> wrote: > > > I've discovered an infelicity (I guess) in qbeta(): it's not a bug, > since there's a clear warning about lack of convergence of the numerical > algorithm ("full precision may not have been achieved"). I can work > around this, but I'm curious why it happens and whether there's a better > workaround -- it doesn't seem to be in a particularly extreme corner of > parameter space. It happens, e.g., for these parameters: > > phi <- 1.1 > i <- 0.01 > t <- 0.001 > shape1 = i/phi ## 0.009090909 > shape2 = (1-i)/phi ## 0.9 > qbeta(t,shape1,shape2) ## 5.562685e-309 > ## brute-force uniroot() version, see below > Qbeta0(t,shape1,shape2) ## 0.9262824 > > The qbeta code is pretty scary to read: the warning "full precision > may not have been achieved" is triggered here: > > https://github.com/wch/r-source/blob/f8d4d7d48051860cc695b99db9be9cf439aee743/src/nmath/qbeta.c#L530 > > Any thoughts? Should I report this on the bug list? > > > A more general illustration: > http://www.math.mcmaster.ca/bolker/misc/qbeta.png > > === > fun <- function(phi,i=0.01,t=0.001, f=qbeta) { > f(t,shape1=i/phi,shape2=(1-i)/phi, lower.tail=FALSE) > } > ## brute-force beta quantile function > Qbeta0 <- function(t,shape1,shape2,lower.tail=FALSE) { > fn <- function(x) {pbeta(x,shape1,shape2,lower.tail=lower.tail)-t} > uniroot(fn,interval=c(0,1))$root > } > Qbeta <- Vectorize(Qbeta0,c("t","shape1","shape2")) > curve(fun,from=1,to=4) > curve(fun(x,f=Qbeta),add=TRUE,col=2) > > ______________________________________________ > R-devel@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-devel -- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Office: A 4.23 Email: pd....@cbs.dk Priv: pda...@gmail.com ______________________________________________ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel