用概率模擬伯努利變量一種b一種b{aover b}使用有偏見的硬幣
誰能告訴我如何模擬, 在哪裡,使用拋硬幣(根據需要多次)?
我正在考慮使用拒絕抽樣,但無法確定。
因為有無數的解決方案,讓我們找到一個有效的解決方案。
這個背後的想法始於實現伯努利變量的標準方法:比較一個統一的隨機變量到參數. 什麼時候, 返回; 否則,返回.
我們可以使用-coin 作為一個統一的隨機數生成器。生成數字均勻地在任何區間內,擲硬幣。當它是頭時,遞歸地生成一個統一的值在第一區間的一部分;當它是尾巴時,遞歸生成從最後區間的一部分。在某些時候,目標間隔會變得非常小,以至於你如何從中選擇一個數字並不重要:這就是遞歸開始的方式。很明顯,這個過程會產生統一的變量(達到任何所需的精度),這很容易通過歸納證明。
這個想法不是有效的,但它導致了一種有效的方法。 因為在每個階段,您都將從某個給定的間隔中抽取一個數字,為什麼不先檢查是否需要繪製呢?如果您的目標值位於此區間之外,則您已經知道隨機值與目標之間的比較結果。 因此,該算法傾向於快速終止。(這可以解釋為問題中要求的拒絕抽樣程序。)
我們可以進一步優化這個算法。 在任何階段,我們實際上都有兩種可以使用的硬幣:通過重新標記我們的硬幣,我們可以使它成為一個有機會的硬幣. 因此,作為預計算,我們可以遞歸地選擇任何重新標記導致終止所需的較低預期翻轉次數的重新標記。 (此計算可能是一個昂貴的步驟。)
例如,將硬幣與模仿伯努利直接變量:平均需要將近十次翻轉。但是如果我們使用硬幣,那麼只需兩次翻轉,我們肯定會完成,預期的翻轉次數只是.
這是詳細信息。
對任何給定的半開區間進行分區進入區間
這定義了兩個轉換和以半開間隔運行。
作為術語,如果是任何一組實數讓表達式
意思是是下界:對所有人. 相似地,方法是一個上限.
寫. (事實上,如果是真實的而不是理性的;我們只要求.)
這是產生變量的算法使用所需的伯努利參數:
- 放和.
- 儘管{拋硬幣生產. 放 增量.}
- 如果然後設置. 否則,設置.
執行
為了說明,這裡是
R
作為函數的算法的實現draw
。它的參數是目標值和間隔, 原來. 它使用輔助功能s
實現. 雖然它不需要,但它也會跟踪拋硬幣的次數。它返回隨機變量、投擲次數和它檢查的最後一個區間。s <- function(x, ab, p) { d <- diff(ab) * p if (x == 1) c(ab[1], ab[1] + d) else c(ab[1] + d, ab[2]) } draw <- function(target, p) { between <- function(z, ab) prod(z - ab) <= 0 ab <- c(0,1) n <- 0 while(between(target, ab)) { n <- n+1; ab <- s(runif(1) < p, ab, p) } return(c(target > ab[2], n, ab)) }
以它的使用和準確性測試為例,以案例為例和. 讓我們畫畫使用算法的值,報告平均值(及其標準誤差),並指示使用的平均翻轉次數。
target <- 0.01 p <- 0.9 set.seed(17) sim <- replicate(1e4, draw(target, p)) (m <- mean(sim[1, ])) # The mean (m - target) / (sd(sim[1, ]) / sqrt(ncol(sim))) # A Z-score to compare to `target` mean(sim[2, ]) # Average number of flips
在這個模擬中翻轉是頭。雖然低於目標, Z 分數不顯著:這種偏差可以歸因於偶然性。平均翻轉次數為——不到十點。如果我們使用了硬幣,平均值是–仍然與目標沒有顯著差異,但僅平均需要翻轉。