Search code examples
rfunctionoptimizationsolver

Securing the interval in uniroot()?


I was wondering if I could make uniroot() work for when q is between nearly 0 and 2?

Currently, uniroot() works only for q > = 3.

P.S. I changed the interval limits, as well as using extendInt = "downX" in uniroot(), none of which worked.

Here is my R code:

f <- function(x, alpha = .05, q = 2, df1 = 3, df2 = 108){
 alpha - suppressWarnings(pf(q, df1, df2, x, lower.tail = FALSE))
}

curve(f)

sapply(c(.05, .95), 
function(i) uniroot(f, interval = c(0, 2e2), alpha = i, q = 2, df1 = 3, df2 = 108)[[1]])

Solution

  • Your problem is that for the values you have specified (alpha=0.05, q=2, df1=3, df2=108), the cumulative distribution function doesn't seem to attain the value you want for any value of the non-centrality parameter (your x).

    f <- function(x, alpha = .05, q = 2, df1 = 3, df2 = 108){
      alpha - suppressWarnings(pf(q, df1, df2, x, 
                                  lower.tail = FALSE))
    }
    

    You can see this most easily by looking at the curve:

    png("c1.png")
    curve(f(x,alpha=0.05,q=2,df1=3,df2=108),
          from=1e-8,to=2e2,log="x")
    dev.off()
    

    enter image description here

    (The curve starts off negative at x=0 and decreases ...)

    Brute force exploration of what happens for varying alpha, ncp, q ...

    library(ggplot2)
    dd <- expand.grid(alpha=seq(0.025,0.975,by=0.025),
                      q=seq(2,4.75,by=0.25),
                      x=emdbook::lseq(1e-2,2e2,length.out=51))
    dd$z <- c(plyr::aaply(dd,1,
                function(z) with(z,f(alpha=alpha,q=q,x=x))))
    ggplot(dd,aes(x,alpha,z=z))+geom_raster(aes(fill=z))+
      facet_wrap(~q,labeller=label_both)+
      scale_x_log10(expand=c(0,0))+
      scale_y_continuous(expand=c(0,0))+
      geom_contour(breaks=0,colour="red")+
      geom_hline(yintercept=c(0.05,0.95),
                 colour="magenta",lty=2)
    ggsave("c2.png")
    

    The red curves are zero contours: the horizontal dashed lines are the transects along which you're trying to find roots. You can see that for low q the transects don't intersect the contour.

    What this means in statistical terms is not something I've thought about ...

    enter image description here