2016-06-10 95 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")開始,但最好的答案絕對是特定於案例的。