最近我們被客戶要求撰寫關(guān)于貝葉斯推斷的研究報(bào)告，包括一些圖形和統(tǒng)計(jì)輸出。

視頻：線性回歸中的貝葉斯推斷與R語(yǔ)言預(yù)測(cè)工人工資數(shù)據(jù)案例

貝葉斯推斷線性回歸與R語(yǔ)言預(yù)測(cè)工人工資數(shù)據(jù)

，時(shí)長(zhǎng)09:58

在這篇文章中，我將擴(kuò)展從數(shù)據(jù)推斷概率的示例，考慮 0 和 1之間的所有（連續(xù)）值，而不是考慮一組離散的候選概率。這意味著我們的先驗(yàn)（和后驗(yàn)）現(xiàn)在是一個(gè)?probability density function?(pdf) 而不是?probability mass function?(pmf)。

我考慮了從數(shù)據(jù)序列推斷 p0，即零的概率：

我使用 p0的不同先驗(yàn)來(lái)解決相同的問(wèn)題，該先驗(yàn)允許 0 和 1?之間的連續(xù)值而不是一組離散的候選值。

概率

我們推理問(wèn)題的起點(diǎn)是?似然?——觀察到的數(shù)據(jù)序列的概率，寫成就像我知道 p0 的值一樣：

?

為了清楚起見(jiàn)，我們可以插入 p0=0.6，并找到給定未知概率值的指定數(shù)據(jù)序列的概率：

概率的更一般形式，不是特定于所考慮的數(shù)據(jù)序列，是

其中 n0 是零的數(shù)量，n1 是考慮的任何數(shù)據(jù)序列 D 中的 1 的數(shù)量。

先驗(yàn) - Beta 分布

我們使用貝塔分布來(lái)表示我們的先驗(yàn)假設(shè)/信息。其數(shù)學(xué)形式是

其中 α0 和 α1 是我們必須設(shè)置的超參數(shù)，反映我們關(guān)于 p0 值的假設(shè)/信息。但是，只需將p0 視為我們想要推斷的參數(shù)——忽略參數(shù)是概率。

請(qǐng)注意?后驗(yàn) pdf 也將是 Beta Distribution，因此值得努力適應(yīng) pdf。

先驗(yàn)均值——?大多數(shù)人想要一個(gè)數(shù)字或點(diǎn)估計(jì)來(lái)表示推理結(jié)果或先驗(yàn)中包含的信息。然而，在貝葉斯推理方法中，先驗(yàn)和后驗(yàn)都是 pdf 或 pmf。獲得點(diǎn)估計(jì)的一種方法是取?相關(guān)參數(shù)相對(duì)于先驗(yàn)或后驗(yàn)的?平均值。例如，對(duì)于 Beta 先驗(yàn)我們得到：

pdf 被歸一化——?這意味著如果我們將 p0 從 0?積分到 1，我們會(huì)得到一個(gè)：

因?yàn)橐韵玛P(guān)系：

就我們而言，最重要的信息是?b?在 0 到 1 區(qū)間上進(jìn)行歸一化，這對(duì)于像 p0 這樣的概率是必要的。

先驗(yàn)假設(shè)和信息可以通過(guò)設(shè)置超參數(shù)來(lái)反映--超參數(shù)α0α0和α1α1影響pdf的形狀，使先驗(yàn)信息的編碼更加靈活。

例如，使用α0=1, α1=1 不能反映p0 的優(yōu)選值。這個(gè)pdf看起來(lái)像

另一個(gè)先驗(yàn)可以指定 α0=5, α1=5，它更接近 p0=1/2 附近的值

最后，我們可以用α0≠α1得到非對(duì)稱的先驗(yàn)，可以看到α0=2，α1=8。

?

關(guān)于設(shè)置超參數(shù)需要記住的一些事情：

• 如果 α0=α1，先驗(yàn)將是對(duì)稱的，先驗(yàn)平均值等于 Eprior[p0]=1/2。
• 如果 α0≠α1，則先驗(yàn)將是不對(duì)稱的，先驗(yàn)平均值不同于 1/21/2。
• 先驗(yàn)的強(qiáng)度與α0+α1的總和有關(guān)。將α總和與數(shù)據(jù)中的n0+n1進(jìn)行比較，將α的視為假數(shù)。這些和的相對(duì)大小控制著先驗(yàn)和似然對(duì)后驗(yàn)形狀的影響。這在下面的Python例子中會(huì)變得很清楚。

累積分布函數(shù) (cdf) ?beta累積分布函數(shù) (cdf)讓我們計(jì)算 p0 小于或等于值 x 的概率。具體來(lái)說(shuō)，cdf定義為：

該積分也被稱為不完全Beta ing積分，并表示為Ix(α0,α1)。

如果我們想知道p0在數(shù)值xl和xh之間的概率，我們可以用cdf來(lái)計(jì)算。
?

不完全 Beta 積分或 cdf 及其逆積分允許從先驗(yàn)或后驗(yàn)計(jì)算置信區(qū)間。使用這些工具，可以說(shuō)p0 的值有 95% 的概率在某個(gè)范圍內(nèi)——同樣，我們將使用 Python 代碼在下面繪制它。

貝塔分布是這個(gè)問(wèn)題的共軛先驗(yàn)--這意味著后驗(yàn)將具有與先驗(yàn)相同的數(shù)學(xué)形式（它也是一個(gè)貝塔分布），并更新了超參數(shù)。這種數(shù)學(xué)上的 "共鳴 "真的很好，讓我們不用MCMC就能做完整的貝葉斯推斷。

現(xiàn)在我們談?wù)勜惾~斯定理和這個(gè)問(wèn)題的后驗(yàn)pdf。

貝葉斯定理和后驗(yàn)

我們的最終目標(biāo)是后驗(yàn)概率密度函數(shù)，結(jié)合似然和先驗(yàn)，在考慮數(shù)據(jù)后對(duì)我們對(duì)p0的知識(shí)做一個(gè)更新的反映。后驗(yàn)pdf的形式是（在這種情況下）。

換句話說(shuō)，這是?給定數(shù)據(jù)序列?D?和先驗(yàn)假設(shè)的?p0?的概率密度，由具有超參數(shù)?(α0,α1)的 Beta pdf 反映。

在這種情況下，貝葉斯定理采用以下形式：

其中后驗(yàn) P(p0|D,α0,α1) 為藍(lán)色，似然 P(D|p0)為黑色，先驗(yàn) P(p0| α0,α1)是紅色的。請(qǐng)注意，歸一化?邊際似然?（上述等式中的分母）現(xiàn)在是一個(gè)積分。

嘗試將貝葉斯定理視為關(guān)于 p0 從?假設(shè)?(α0,α1) 更新到?假設(shè) + 數(shù)據(jù)?(D,α0,α1) 的信息：

試著把貝葉斯定理看作是關(guān)于p0的信息被從假設(shè)（α0,α1）更新為假設(shè)+數(shù)據(jù)（D,α0,α1）。

為了得到后驗(yàn)pdf，我們必須在貝葉斯定理的分母上做積分。在這種情況下，利用貝塔分布的特性，就可以進(jìn)行計(jì)算。該積分如下。

最后一行的積分定義了一個(gè)貝塔函數(shù)，在關(guān)于先驗(yàn)的一節(jié)中討論過(guò)，并且有一個(gè)已知的結(jié)果。
?

這意味著分母，也叫邊際似然，等于。?

同樣，我們得到這個(gè)結(jié)果是因?yàn)锽eta分布是我們所考慮的伯努利過(guò)程概率的共軛先驗(yàn)。請(qǐng)注意，來(lái)自先驗(yàn)的超參數(shù)已經(jīng)被計(jì)數(shù)數(shù)據(jù)所更新。

這與人們預(yù)期的完全一樣，不需要做所有的數(shù)學(xué)計(jì)算。在任何情況下，在用Python實(shí)現(xiàn)這一點(diǎn)之前，有幾個(gè)注意事項(xiàng)。
?

• 后驗(yàn)pdf在0到1的區(qū)間內(nèi)被歸一化，就像我們推斷p0這樣的概率時(shí)需要的那樣。
• 后驗(yàn)平均數(shù)，是對(duì)我們的推斷給出一個(gè)點(diǎn)估計(jì)的方法是

• 后驗(yàn)的cdf和先驗(yàn)的一樣，因?yàn)槲覀內(nèi)匀挥幸粋€(gè)Beta分布--只不過(guò)，現(xiàn)在的參數(shù)是用數(shù)據(jù)更新的。在任何情況下，我們都可以用不完整的Beta積分和它的逆向找到置信區(qū)間，如上所述。

Python 中的推理代碼

首先，我們導(dǎo)入一些包，使用這些包來(lái)計(jì)算和繪制先驗(yàn)、似然和后驗(yàn)。此外，使用 matplotlib，在本例中為 ggplot，創(chuàng)建漂亮的圖。

概率

def _int(sef, daa):""二進(jìn)制數(shù)據(jù)的概率。""elfus = {s:0 for s in ['0', '1']}def proes_data(sef, ata):""處理數(shù)據(jù)。""for s in ['0', '1']:unts[s] =cont(s)""處理概率。""n0 = couts['0']n1 = couts['1']if p0 != 0 and p0 != 1:# 典型情況lgpdaa = n0*log(p0) + \n1*log(1.-p0)rdta = exp(lgpata)elif p0 == 0 and n0 != 0:# 如果n0不是0，p0就不可能是0lpr_dta = -np.infp_dta = exp(lgrata)elif p0 == 0 and n0 == 0:# 數(shù)據(jù)與p0=0一致lgpr_at = n1*np.log(1.-p0)p_ata =.exp(lo_dta)elif p0 == 1 and n1 != 0:# 如果n1不是0，p0就不可能是1lpdta = -infpraa = exp(lor_ta)elif p0 == 1 and n1 == 0:# 數(shù)據(jù)與p0=1一致log_dta = n0*log(p0)def prob(sef, p0):""獲取數(shù)據(jù)的概率。"""prdat,_ = ._procs_proabiti(p0)def logprob(sef, p0):""獲取數(shù)據(jù)的概率對(duì)數(shù)。"""_, lgpr_a = sel_presplies(p0)

先驗(yàn)分布

我們的先驗(yàn)類基本上是一個(gè)圍繞 scipy的包，有一個(gè)繪圖方法。注意 plot() 方法得到了 Beta 分布的平均值，并使用? interval() 方法得到了一個(gè)概率為 95% 的區(qū)域--這是使用不完整的 Beta 積分和上面討論的它的逆值完成的。
?

def _int_(sef, alpa0=1, alph1=1):""二進(jìn)制數(shù)據(jù)的β先驗(yàn)。""elfa0 = alha0sl.1 = alha1slf.0r = ba(self.0, sef.a1)def intrval(slf, pb):""包含`prob'的pdf區(qū)域的端點(diǎn)。例如：interval(0.95)"""return sef.pvitervl(rob)def mean(self):"""返回先驗(yàn)平均數(shù)。"""def pdf(self, p0):"""p0處的概率密度。"""return self.p.pdf(p0)def plot(self):""顯示平均值和95%置信區(qū)間的圖。""fig, ax = tsuplots(1, 1)x = np.rae(0., 1., 0.01)# 獲得先驗(yàn)平均數(shù) p0mean = mean()# 獲得包含95%概率的低/高點(diǎn)low_p0, high_p0 = l.nterval(.95)xo =nang(_0, hih_p0, 0.01)# 繪制pdfax.plot(x, self.pdf(x), 'r-')# 填充95%的區(qū)域between(x_prob 0, sel.pdf(x_pob)# 平均值ax.stem([mean], [sf.df(mea)

讓我們使用新代碼繪制一些具有一序列參數(shù)的 Beta pdf。

統(tǒng)一先驗(yàn)

帶點(diǎn)的垂直線顯示 pdf 均值的位置。陰影區(qū)域表示對(duì)于給定的 α0 和 α1 值，概率為 95% 的（對(duì)稱）區(qū)域。如果您想要平均值和置信區(qū)間的實(shí)際值，也可以獲取這些值：

print("先驗(yàn)均值: {}".format(pri.mean()))

上面的其他先前示例也有效：

prior(5, 5) plot()

和

por(2, 8) iplot()

?

了解超參數(shù)所反映的先前假設(shè)的均值和不確定性很有用。

后驗(yàn)

最后，我們?yōu)楹篁?yàn)構(gòu)建類。正如您所料，我將數(shù)據(jù)和先驗(yàn)作為參數(shù)，并從這些元素中提取后驗(yàn)所需的參數(shù)。

def __int_(sf datprior):""后驗(yàn)。data: 一個(gè)數(shù)據(jù)樣本的列表先驗(yàn)：β先驗(yàn)類的一個(gè)實(shí)例"""elflilihod = liliod(dta)sefprir = rioself._css_steror()def _pces_posteror(slf):""使用傳遞的數(shù)據(jù)和先驗(yàn)來(lái)處理后驗(yàn)。""# 從似然和先驗(yàn)中提取n0, n1, a0, a1seln = slfliklihod.counts['0']sel.n1 = elf.lkelihodcnts['1']lfa0 = sf.prir.0self.a = sef.priorael0rv= beta(selfa0 + sfn0,sef.a1 + slf.1)def interalself, prob):""含`prob`的pdf區(qū)域的端點(diǎn)。例如：interval(0.95)"""def mean(sef):"""返回后驗(yàn)平均數(shù)。"""def pdf(sef, p0):"""p0處的概率密度。"""def plot(slf):""顯示先驗(yàn)、似然和后驗(yàn)的圖。""## 先驗(yàn)# 得到先驗(yàn)平均數(shù) p0pri_mean =eorman()#得到包含95%概率的低/高分值plo_p0,= interval(0.95)prob = arange(low_p0, i_hgh_0, 0.01)# 繪制pdfplot(prior.pdf(x)# 填充95%的區(qū)域x.ll_between(pri_p )# 平均值astm([pri_mean])## 似然# 繪制似然圖li = [sel.likliood.pro]# ##后驗(yàn)#獲得后驗(yàn)平均數(shù)p0ostmen = mean()#得到包含95%概率的低/高點(diǎn)ow_p0, _high_p0 = interval(0.95)prob = np.rngest_low_p0po_highp0 0.01)# 繪制pdfplot(x, slf.pd(x)# 填充95%的區(qū)域fil_etween(pos_xob, 0,self.df(pt_pr )# 平均值ax2].t_mean])

基本代碼就是這樣，讓我們??做一些例子。

例子

讓我們從數(shù)據(jù)和統(tǒng)一先驗(yàn)的示例開始。

# 數(shù)據(jù) data1 = [0,0,0,0,1,1,0,0,0,1]# 先驗(yàn) pro1 = prir(1, 1)# 后驗(yàn) ot1 = psteior(daa1, ior1) plot()

這里需要注意的事項(xiàng)：?

• 先驗(yàn)是統(tǒng)一的。這意味著概率和后驗(yàn)具有相同的形狀。
• 95%的置信區(qū)間同時(shí)顯示在先驗(yàn)和后驗(yàn)中。

接下來(lái)，讓我們考慮具有不統(tǒng)一先驗(yàn)的相同數(shù)據(jù)。數(shù)據(jù)集長(zhǎng)度為 10，因此 n0+n1=10。讓我們用 α0+α1=10 ，設(shè)置先驗(yàn)，但先驗(yàn)在與似然不同的位置達(dá)到峰值（也許有專家說(shuō)這應(yīng)該是先驗(yàn)設(shè)置）：

# 先驗(yàn) prir(4, 6)# 后驗(yàn) ps2 = postior(da1, pio2) plot()

顯然數(shù)據(jù)和專家在這一點(diǎn)上存在分歧。然而，因?yàn)橄闰?yàn)的權(quán)重設(shè)置為 10 并且數(shù)據(jù)序列的長(zhǎng)度為 10，所以后驗(yàn)峰值位于先驗(yàn)峰值和似然峰值的中間。嘗試使用這種效果來(lái)更好地理解先驗(yàn)超參數(shù)、數(shù)據(jù)集長(zhǎng)度和結(jié)果后驗(yàn)之間的相互作用。

作為最后一個(gè)例子，我們考慮最后一個(gè)例子的兩個(gè)變體，?首先我們使用統(tǒng)一先驗(yàn)：

# 設(shè)定概率為0 p0 = 0.2 #設(shè)置rng種子為2 n.rand.sed(2) # 產(chǎn)生數(shù)據(jù) dta2 =rando.i([0,1], 500, p=[p0, 1.-p0])# 先驗(yàn) prior(1,1)# 后驗(yàn) poteir(daa2, pio3)

請(qǐng)注意，概率和后驗(yàn)峰值在同一個(gè)地方，正如我們所期望的那樣。但是，由于數(shù)據(jù)集較長(zhǎng)（500 個(gè)值），峰值要強(qiáng)得多。

最后，我們?cè)谕粩?shù)據(jù)集上使用“錯(cuò)誤先驗(yàn)”。在這種情況下，我們將保持先驗(yàn)強(qiáng)度為 10，即 α0+α1=10：

# 先驗(yàn) prior(6,4) # 后驗(yàn) poseor(data, pior4)

請(qǐng)注意，盡管先驗(yàn)在錯(cuò)誤的位置達(dá)到峰值，但概率和后驗(yàn)非常相似。這個(gè)例子表明，如果先驗(yàn)沒(méi)有設(shè)置得太強(qiáng)，合理數(shù)量的數(shù)據(jù)應(yīng)該產(chǎn)生不錯(cuò)的推理。一般來(lái)說(shuō)，最好讓 n0+n1>α0+α1 并考慮先驗(yàn)和后驗(yàn)的形狀。

---

總結(jié)

以上是生活随笔為你收集整理的PYTHON贝叶斯推断计算：用BETA先验分布推断概率和可视化案例的全部?jī)?nèi)容，希望文章能夠幫你解決所遇到的問(wèn)題。

如果覺(jué)得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。