2016-06-10 54 views
2

我试图用optimize()找n个以下函数的最小值(Clopper皮尔森下限):优化全球最小

f <- function (n, p=0.5) 
(1 + (n - p*n + 1)/
    (p*n*qf(p= .025, df1= 2*p, df2= 2*(n - p + 1))))^-1 

这里就是我试图优化它:

n_clop <- optimize(f.1, c(300,400), maximum = FALSE, p=0.5) 
n_clop 

我在间隔[300,400]上做了这个,因为我怀疑它的值在它之内,但最终我想在0和无穷之间进行优化。看起来这个命令产生了一个局部最小值,因为不管它产生的时间间隔的下限是否为最小值 - 这不是我从clopper-pearson得知的。所以,我的两个问题是如何正确地找到R中的全局最小值以及如何在任何时间间隔内如此超过?

+0

此外,我希望函数等于0.5(比例置信区间的半宽)并计算n。不确定我的设置是否正确。 –

+0

顺便说一下,SO政策是你不需要在问题标题中包含编程语言(R)的名称 - [r]标记应该足够了 –

+0

没有通用的方法来执行此操作每个功能。你愿意承担什么来解决这个问题?全局优化对任何语言都不是一个容易解决的问题。 – MrFlick

回答

2

很简单望了Wikipedia page you linked并没有看到您的公式中的任何明显的错别字(虽然我觉得它应该是0.975 = 1-α/ 2,而不是0.025 =阿尔法/ 2?) 。然而,评估你在很大范围内编码的功能意味着没有局部最小值让你感到困惑。我强烈的猜测是,或者你的逻辑错误(即n - > 0的确是正确的答案),或者你没有编码你认为自己编码的内容,这是由于一个错字(可能在维基百科文章中,尽管这似乎不太可能)或一个thinko。

f <- function (n, p=0.5) 
(1 + (n - p*n + 1)/
    (p*n*qf(p= .025, df1= 2*p, df2= 2*(n - p + 1))))^-1 

确认你得到正确的答案选择的时间间隔:

curve(f(x),c(300,400)) 

评估在广泛的范围内(N = 0.00001到1000000):

curve(f(10^x),c(-5,7)) 

enter image description here

正如@MrFlick所言,全局优化很难。您可以从optim(...method="SANN")开始,但最好的答案绝对是特定于案例的。