2016-11-20 81 views
2

有些情況下,我並不真正對貝葉斯推斷的完整後驗感興趣,而只是最大可能性(或對於適當選擇的先驗的最大後驗),並且可能是Hessian。 PyMC3有這樣做的功能,但find_MAP似乎返回轉換形式的模型參數,取決於它們的先前分佈。有沒有一種簡單的方法來從這些未轉換的值? find_hessian的輸出對我來說更不明顯,但它最有可能也在轉換的空間中。使用PyMC3計算最大似然率

回答

1

可能是簡單的解決方案將是傳遞參數transform=None,避免PyMC3執行轉換,然後使用find_MAP

我讓你,比如一個簡單的模型。

data = np.repeat((0, 1), (3, 6)) 
with pm.Model() as normal_aproximation: 
    p = pm.Uniform('p', 0, 1, transform=None) 
    w = pm.Binomial('w', n=len(data), p=p, observed=data.sum()) 
    mean_q = pm.find_MAP() 
    std_q = ((1/pm.find_hessian(mean_q))**0.5)[0] 
print(mean_q['p'], std_q) 

您是否考慮過使用ADVI

+0

我有,但在我的情況下,模型很簡單,MAP + hessian很好地描述了後驗分佈(它接近於正常)。 MAP + hessian比ADVI或MC-sampling要快得多。當然,這是濫用pymc3(通過不採樣),但他們的模型規範是偉大的,我想保留採樣的選項,以防模型/數據的變化導致更復雜的後驗。 – burnpanck

0

我再次遇到這個問題,並找到了從轉換後的值中獲取未轉換值的方法。以防萬一別人需要這樣做。它的要點是未轉換的值本質上是可以在給定轉換後的值的情況下評估的表達式。通過提供Model.fn()函數,PyMC3在這方面有所幫助,該函數創建了一個按名稱接受值的評估函數。現在,您只需將未轉換的感興趣變量提供給outs參數。一個完整的例子:

data = np.repeat((0, 1), (3, 6)) 
with pm.Model() as normal_aproximation: 
    p = pm.Uniform('p', 0, 1) 
    w = pm.Binomial('w', n=len(data), p=p, observed=data.sum()) 
    map_estimate = pm.find_MAP() 
# create a function that evaluates p, given the transformed values 
evalfun = normal_aproximation.fn(outs=p) 
# create name:value mappings for the free variables (e.g. transformed values) 
inp = {v:map_estimate[v.name] for v in model.free_RVs} 
# now use that input mapping to evaluate p 
p_estimate = evalfun(inp) 

outs還可以接收的變量的列表,evalfun然後將輸出以相同的順序對應的變量的值。