给定R中的条件，如何向后求解未知数？

Question

对于一组输入（见下文），当我运行以下 R 代码时，我得到两个存储在 ncp 中的答案。但我想知道df2 应该是什么，以便ncp（即abs(ncp[2] - ncp[1])）中这两个答案之间的区别是.15？

所以，其他一切都已解决，df2 应该是什么，以便abs(ncp[2] - ncp[1]) = .15？这可以在 R 中完成吗？

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))
}
# Example of use:
 ncp(df2 = 108) # Answers: 0.1498627 0.4100823
                # What should `df2` be so that the difference between 2 answers is `.15`

Answer 1

stats 包中的optim 和optimize 函数是处理此类问题的好方法。 optimize 更简单，处理单变量情况。这是一个示例解决方案：

# Define a function that we want to minimise the output of
ncp_diff <- function(df2, target = 0.15){
  the_ncp <- ncp(df2)
  return(abs(abs(the_ncp[2] - the_ncp[1]) - target))
}

# try it out - how far out is df2 = 108 from our target value?
ncp_diff(108) # 0.1102

# find the value of ncp_diff's argument that minimises its return:
optimize(ncp_diff, interval = c(0, 1000)) # minimum = 336.3956

ncp(336.3956) # 0.218, 0.368

请注意，虽然 336.3956 是 a 解决方案，但它不一定是 解决方案。您需要注意局部最小值。不过，处理这个问题超出了简单答案的范围。