具有numpy的泊松置信區間

Question

我想將泊松連續誤差條放在我用matplotlib製作的直方圖上，但我似乎無法找到一個numpy函數，它會給出95％的置信區間，假設泊松數據。理想情況下，解決方案不依賴於scipy，但任何東西都可以工作。這樣的功能是否存在？我已經發現了很多關於bootstrapping的內容，但在我看來這似乎有點過分。

Answer 1

使用scipy.stats.poisson和interval方法：

>>> scipy.stats.poisson.interval(0.95, [10, 20, 30]) 
(array([ 4., 12., 20.]), array([ 17., 29., 41.])) 

即使它不僅使有限的意義上，以計算用於非整數值的泊松分佈，由OP要求確切置信區間可以被計算這是可以做到如下：

>>> data = np.array([10, 20, 30]) 
>>> scipy.stats.poisson.interval(0.95, data) 
(array([ 4., 12., 20.]), array([ 17., 29., 41.])) 
>>> np.array(scipy.stats.chi2.interval(.95, 2 * data))/2 - 1 
array([[ 3.7953887 , 11.21651959, 19.24087402], 
     [ 16.08480345, 28.67085357, 40.64883744]]) 

---

它也可以使用ppf方法：

>>> data = np.array([10, 20, 30]) 
>>> scipy.stats.poisson.ppf([0.025, 0.975], data[:, None]) 
array([[ 4., 17.], 
     [ 12., 29.], 
     [ 20., 41.]]) 

但由於分佈是離散的返回值將是整數，置信區間不會跨越正好95％：

>>> scipy.stats.poisson.ppf([0.025, 0.975], 10) 
array([ 4., 17.]) 
>>> scipy.stats.poisson.cdf([4, 17], 10) 
array([ 0.02925269, 0.98572239]) 

Answer 2

我最後寫基於我自己的函數some properties I found on Wikipedia。

def poisson_interval(k, alpha=0.05): 
    """ 
    uses chisquared info to get the poisson interval. Uses scipy.stats 
    (imports in function). 
    """ 
    from scipy.stats import chi2 
    a = alpha 
    low, high = (chi2.ppf(a/2, 2*k)/2, chi2.ppf(1-a/2, 2*k + 2)/2) 
    if k == 0: 
     low = 0.0 
    return low, high 

這將返回連續的（而不是離散的）邊界，這在我的字段中更爲標準。

Answer 3

此問題出現在天文學很多（我的領域！），這紙是去到這些置信區間參考：Gehrels 1980

它有很多數學在它的任意置信區間與泊松統計量，但對於雙側95％置信區間（對應於2西格瑪高斯置信區間，或本文中上下文中的S = 2），當N個事件發生時，上下置信區間的一些簡單分析公式被測量的是

upper = N + 2. * np.sqrt(N + 1) + 4./3. 
lower = N * (1. - 1./(9. * N) - 2./(3. * np.sqrt(N))) ** 3. 

我在哪裏把它們放在Python格式中爲你alr伊迪。所有你需要的是numpy或你最喜歡的平方根模塊。請記住，這些會給你事件的上限和下限 - 而不是+/-值。你只需從這兩個中減去N就可以得到這些。

有關這些公式的準確性，請查閱您需要的置信區間，但這些公式對於大多數實際應用而言應該足夠精確。