CODEWITHSUNDEEP
rfunctionoptimization
rnorouzian
rnorouzian

Reputation: 7517

optimize() giving different answers?

Below, I'm trying to solve for ncp (there is one answer). But I'm wondering why when I extend the interval argument in optimize the answer drastically changes?

Could I use uniroot instead of optimize here?

f <- function(pwr, q, df1, df2, ncp){
 abs(pwr - pf(q, df1, df2, ncp, lower.tail = FALSE))
}

optimize(f, interval = c(0, 1e2), pwr = .8, q = 2.5, df1 = 3, df2 = 108)[[1]]
# [1] 10.54639  !!! HERE

optimize(f, interval = c(0, 5e2), pwr = .8, q = 2.5, df1 = 3, df2 = 108)[[1]]
# [1] 499.9999  !!! HERE

Upvotes: 1

Views: 61

Answers (1)

Ben Bolker
Ben Bolker

Reputation: 226057

Because the rightmost part of the curve is too flat - all values beyond 150 are identical.

Utility function:

f2 <- function(x) f(x, pwr = .8, q = 2.5, df1 = 3, df2 = 108)
cc <- curve(f2(x)-0.2,from=150,to=500)
unique(cc$y)
## [1] -5.551115e-17

uniroot() does indeed work fine: we have to change the function f to return a signed value .

f <- function(pwr, q, df1, df2, ncp){
   pwr - pf(q, df1, df2, ncp, lower.tail = FALSE)
}
uniroot(f, interval = c(0, 5e2), pwr = .8, q = 2.5, df1 = 3, df2 = 108)
## $root
## [1] 10.54641
## $f.root
## [1] -3.806001e-08 
## etc.

In general, converting root-finding problems to minimum-finding problems by squaring or taking the absolute value is a fragile strategy (I read about this in Numerical Recipes years ago ...)

Upvotes: 3

Related Questions