我们可以在 R 中使用 uniroot() 而不是 optimize() 来进行这种最小化吗？

Question

在我下面的 R 代码的最后一行中，我使用 optimize() 找到最小化 ncp_diff 函数的 df2。

但是，我想知道我是否可以 uniroot() 而不是 optimize() 来进行这种最小化？

alpha = c(.025, .975); df1 = 3; peta = .3   # The input

f <- function(alpha, q, df1, df2, ncp){     # Notice `ncp` is the unknown
  alpha - suppressWarnings(pf(q = (peta / df1) / ((1 - peta)/df2), df1, df2, ncp, lower = FALSE))
}

ncp <- function(df2){      # Root finding: finds 2 `ncp` for a given `df2`

 b <- sapply(c(alpha[1], alpha[2]),
       function(x) uniroot(f, c(0, 1e7), alpha = x, q = peta, df1 = df1, df2 = df2)[[1]])

 b / (b + (df2 + 4))
}

ncp_diff <- function(df2, target = 0.15){
 the_ncp <- ncp(df2)
  return(abs(abs(the_ncp[2] - the_ncp[1]) - target))
 }

optimize(ncp_diff, c(0, 1000)) ## HERE can I use `uniroot()` instead of `optimize()`

Answer 1

alpha = c(.025, .975); df1 = 3; peta = .3   # The input

f <- function(alpha, q, df1, df2, ncp){     # Notice `ncp` is the unknown
    alpha - suppressWarnings(pf(q = (peta / df1) / ((1 - peta)/df2), df1, df2, ncp, lower = FALSE))
}

ncp <- function(df2){      # Root finding: finds 2 `ncp` for a given `df2`

    b <- sapply(c(alpha[1], alpha[2]),
                function(x) uniroot(f, c(0, 1e7), alpha = x, q = peta, df1 = df1, df2 = df2)[[1]])

    b / (b + (df2 + 4))
}

ncp_diff <- function(df2, target = 0.15){
    the_ncp <- ncp(df2)
    return((the_ncp[2] - the_ncp[1]) - target)
}

uniroot(ncp_diff, c(100, 1000)) #
$root
[1] 336.3956

$f.root
[1] 3.74663e-09

$iter
[1] 7

$init.it
[1] NA

$estim.prec
[1] 6.103516e-05

编辑：

为了使用相同的间隔 (0,1000)，我们可以寻找一种方法来解决下限值和上限值都在数字线的同一侧产生结果的情况。由于这是 r 中的一个错误，我们可以通过 tryCatch

ncp <- function(df2){      # Root finding: finds 2 `ncp` for a given `df2`

  b <- sapply(c(alpha[1], alpha[2]),
              function(x) 
    tryCatch(uniroot(f, c(0, 1e7), alpha = x, q = peta, df1 = df1, df2 = df2)[[1]],
             error =function(e)NA ))
  if(any(is.na(b)))b= c(1,10)

  b / (b + (df2 + 4))
}


uniroot(ncp_diff, c(0, 1000)) #
$root
[1] 336.3956

$f.root
[1] -2.132438e-09

$iter
[1] 8

$init.it
[1] NA

$estim.prec
[1] 6.103516e-05

我们可以在 R 中使用 uniroot() 而不是 optimize() 来进行这种最小化吗？

Could we use uniroot() instead of optimize() for this minimization in R?

math

optimization

r

function

linear-algebra