python實現貝葉斯推斷的例子

 更新時間:2021年09月10日 11:31:03   作者:笨牛慢耕  
本文介紹一個貝葉斯推斷的python實現,并展現了基于標量運算的實現和基于numpy的矩陣運算的實現之間的差別,感興趣的可以了解一下

1. 前言

        本文介紹一個貝葉斯推斷的python實現例,并展現了基于標量運算的實現和基于numpy的矩陣運算的實現之間的差別。

2. 問題描述

本問題例取自于Ref1-Chapter1.

問題描述:假設有一個制作燈泡的機器。你想知道機器是正常工作還是有問題。為了得到答案你可以測試每一個燈泡,但是燈泡數量很多,每一個都測試在實際生產過程中可能是無法承受的。使用貝葉斯推斷,你可以基于少量樣本(比如說抽檢結果)來估計機器是否在正常地工作(in probabilistic way)。

構建貝葉斯推斷時,首先需要兩個要素:

(1) 先驗分布

(2) 似然率

先驗分布是我們關于機器工作狀態的初始信念。首先我們確定第一個刻畫機器工作狀態的隨機變量,記為M。這個隨機變量有兩個工作狀態:{working, broken},以下簡寫成{w, br}(縮寫成br是為了與下面的Bad縮寫成b區分開來)。作為初始信念,我們相信機器是好的,是可以正常工作的,定義先驗分布如下:

P(M=working) = 0.99

P(M=broken ) = 0.01

這表明我們對于機器正常工作的信念度很高,有99%的概率能夠正常工作。

第二個隨機變量是L,表示機器生產的燈泡的工作狀態。燈泡可能是好,也可能是壞的,包含兩個狀態:{good, bad},以下簡寫成{g,b},注意br與b的區別。

我們需要基于機器工作狀態給出L的先驗分布,也就是條件概率P(L|M),在貝葉斯公式中它代表似然概率(likelihood)。

定義這個似然概率分布(由于M和L各有兩種狀態,所以一共包含4個條件概率)如下:

P(L=Good|M=w) = 0.99

P(L=Bad |M=w) = 0.01

P(L=Good|M=br ) = 0.6

P(L=Bad |M=br ) = 0.4

以上似然概率表明,在機器正常時我們相信每生成100個燈泡只會有一個壞的,而機器不正常時也不是所有燈泡都是壞的,而是有40%會是壞的。為了實現的方便,可以寫成如下的矩陣形式:

        現在,我們已經完整地刻畫了貝葉斯模型,可以用它來做一些神奇的估計和預測的工作了。

        我們的輸入是一些燈泡的抽檢結果。假設我們抽檢了十個燈泡其抽檢結果如下:

   {bad, good, good, good, good, good, good, good, good, good}

        讓我們來看看基于貝葉斯推斷的我們對于機器工作狀態的信念(后驗概率)如何變化。

3. 貝葉斯規則

        貝葉斯推斷規則以貝葉斯公式的形式表示為:

         具體映射到本問題中可以表達如下:

         貝

葉斯推斷的優先在于可以以在線(online)的方式進行,即觀測數據可以一個一個地到來,每次受到一個新的觀測數據,就進行一次基于貝葉斯公式的后驗概率的計算更新,而更新后的后驗概率又作為下一貝葉斯推斷的先驗概率使用。因此在線的貝葉斯推斷的基本處理流程如下所示:

4. Bayes engine: scalar implementation

        首先,我們以標量運算的方式寫一個函數來進行bayes推斷處理。

        prior以向量的形式存儲先驗概率分布,prior[0]表示P(M=working),prior[1]表示P(M=broken)。

        likelihood以矩陣的形式方式存儲似然概率分布。其中第1行表示P(L/M=working),第2行表示P(L/M=broken).

        在本例中,當輸入 時,evidence的計算式(注意evidence是依賴于輸入的觀測數據的)是:

        注意,當我寫P(w)其實是表示P(M=w),而 其實是表示 ,余者類推。根據上下文,這些應該不會導致混淆。

        第一個函數的代碼如下:

def bayes_scalar(prior, likelihood, data):
    """
    Bayesian inference function example.
    Parameters
    ----------
    prior : float, 1-D vector
        prior information, P(X).
    likelihood : float 2-D matrix
        likelihood function, P(Y|X).
    data : List of strings. Value: 'Good','Bad'
        Observed data samples sequence
    Returns
    -------
    posterior : float
        P(X,Y), posterior sequence.
    """
    
    posterior = np.zeros((len(data)+1,2))
    posterior[0,:] = prior  # Not used in computation, just for the later plotting
 
    for k,L in enumerate(data):
        if L == 'good':
            L_value = 0
        else:
            L_value = 1
        #print(L, L_value, likelihood[:,L_value])
 
        evidence      = likelihood[0,L_value] * prior[0] + likelihood[1,L_value] * prior[1]
        LL0_prior_prod= likelihood[0,L_value] * prior[0]                
        posterior[k+1,0]  = LL0_prior_prod / evidence
 
        LL1_prior_prod= likelihood[1,L_value] * prior[1]        
        posterior[k+1,1]  = LL1_prior_prod / evidence
        
        prior = posterior[k+1,:] # Using the calculated posterior at this step as the prior for the next step
                
    return posterior

 5. Bayes engine: vectorization

        我們注意到,evidence的計算可以表示成兩個向量的點積,如下所示。

        這樣就非常方便用numpy來實現了。本例中每個隨機變量只有兩種取值,在復雜的情況下,每個隨機變量有很多種取值時,有效利用向量或矩陣的運算是簡潔的運算實現的必不可缺的要素。以上這兩個向量的點積可以用numpy.dot()來實現。

        另外,likelihood和prior的乘積是分別針對M的兩種狀態進行計算(注意,我們需要針對M的兩種不同狀態分別計算posterior),不是用向量的點積進行計算,而是一種element-wise multiplication,可以用numpy.multiply()進行計算。所以在vectorization版本中貝葉斯更新處理削減為兩條語句,與上面的scalar版本相比顯得非常優雅簡潔(好吧,也許這個簡單例子中還顯不出那么明顯的優勢,但是隨著問題的復雜度的增加,這種優勢就會越來越明顯了。)

         由此我們得到向量化處理的函數如下: 

def bayes_vector(prior, likelihood, data):
    """
    Bayesian inference function example.
    Parameters
    ----------
    prior : float, 1-D vector
        prior information, P(X).
    likelihood : float 2-D matrix
        likelihood function, P(Y|X).
    data : List of strings. Value: 'Good','Bad'
        Observed data samples sequence
    Returns
    -------
    posterior : float
        P(X,Y), posterior sequence.
    """
    
    posterior = np.zeros((len(data)+1,2))
    posterior[0,:] = prior  # Not used in computation, just for the later plotting
 
    for k,L in enumerate(data):
        if L == 'good':
            L_value = 0
        else:
            L_value = 1
        #print(L, L_value, likelihood[:,L_value])
 
        evidence          = np.dot(likelihood[:,L_value], prior[:])
        posterior[k+1,:]  = np.multiply(likelihood[:,L_value],prior)/evidence
 
        prior = posterior[k+1,:] # Using the calculated posterior at this step as the prior for the next step
        
    return posterior

6. 測試

        讓我們來看看利用以上函數對我們的觀測數據進行處理,后驗概率將會如何變化。

import numpy as np
import matplotlib.pyplot as plt
 
prior      = np.array([0.99,0.01])
likelihood = np.array([[0.99,0.01],[0.6,0.4]])
data       = ['bad','good','good','good','good','good','good','good']        
 
posterior1 = bayes_scalar(prior,likelihood,data)
posterior2 = bayes_vector(prior,likelihood,data)
 
if np.allclose(posterior1,posterior2):
    print('posterior1 and posterior2 are identical!')
 
fig, ax = plt.subplots()
ax.plot(posterior1[:,0])
ax.plot(posterior1[:,1])
ax.grid()
# fig.suptitle('Poeterior curve vs observed data')
ax.set_title('Posterior curve vs observed data')
plt.show()

        運行以上代碼可以得到后驗概率的變化如下圖所示(注意第一個點是prior):

         當然以上代碼也順便驗證了一下兩個版本的bayes函數是完全等價的。

7. 后記

        大功告成。

        第1個關于貝葉斯統計的學習的程序和第1篇關于貝葉斯統計的學習的博客。

        其它有的沒的等想到了什么再回頭來寫。

[Ref1] 《概率圖模型:基于R語言》,David Bellot著, 魏博譯

到此這篇關于python實現貝葉斯推斷的例子的文章就介紹到這了,更多相關python 貝葉斯推斷內容請搜索腳本之家以前的文章或繼續瀏覽下面的相關文章希望大家以后多多支持腳本之家!

相關文章

  • python xlwt如何設置單元格的自定義背景顏色

    python xlwt如何設置單元格的自定義背景顏色

    這篇文章主要介紹了python xlwt如何設置單元格的自定義背景顏色,文中通過示例代碼介紹的非常詳細,對大家的學習或者工作具有一定的參考學習價值,需要的朋友可以參考下
    2019-09-09
  • 淺析Python實現DFA算法

    淺析Python實現DFA算法

    DFA全稱為Deterministic Finite Automaton,即確定有窮自動機。特征:有一個有限狀態集合和一些從一個狀態通向另一個狀態的邊,每條邊標記有一個符號,其中一個狀態是初態,某些狀態是終態。不同于不確定的有限自動機,DFA中不會有從同一狀態出發的兩條邊標志有相同的符號
    2021-06-06
  • Django項目使用CircleCI的方法示例

    Django項目使用CircleCI的方法示例

    這篇文章主要介紹了Django項目使用CircleCI的方法示例,文中通過示例代碼介紹的非常詳細,對大家的學習或者工作具有一定的參考學習價值,需要的朋友們下面隨著小編來一起學習學習吧
    2019-07-07
  • 使用Python實現BT種子和磁力鏈接的相互轉換

    使用Python實現BT種子和磁力鏈接的相互轉換

    這篇文章主要介紹了使用Python實現BT種子和磁力鏈接的相互轉換的方法,有時比如迅雷無法加載磁力鏈接或者無法上傳附件分享時可以用到,需要的朋友可以參考下
    2015-11-11
  • Python基于select實現的socket服務器

    Python基于select實現的socket服務器

    這篇文章主要介紹了Python基于select實現的socket服務器,實例分析了Python基于select與socket模塊實現socket通信的相關技巧,需要的朋友可以參考下
    2016-04-04
  • Python collections模塊的使用技巧

    Python collections模塊的使用技巧

    Python的最大優勢之一是其廣泛的模塊和軟件包。這將Python的功能擴展到許多受歡迎的領域,包括機器學習、數據科學和Web開發等, 其中最好的模塊之一是Python的內置collections 模塊。
    2021-04-04
  • python seaborn heatmap可視化相關性矩陣實例

    python seaborn heatmap可視化相關性矩陣實例

    這篇文章主要介紹了python seaborn heatmap可視化相關性矩陣實例,具有很好的參考價值,希望對大家有所幫助。一起跟隨小編過來看看吧
    2020-06-06
  • python實現四人制撲克牌游戲

    python實現四人制撲克牌游戲

    這篇文章主要介紹了python實現四人制撲克牌游戲,文中示例代碼介紹的非常詳細,具有一定的參考價值,感興趣的小伙伴們可以參考一下
    2020-04-04
  • Python的垃圾回收機制詳解

    Python的垃圾回收機制詳解

    這篇文章主要介紹了Python的垃圾回收機制詳解,文中通過示例代碼介紹的非常詳細,對大家的學習或者工作具有一定的參考學習價值,需要的朋友可以參考下
    2019-08-08
  • 通過Python pyecharts輸出保存圖片代碼實例

    通過Python pyecharts輸出保存圖片代碼實例

    這篇文章主要介紹了通過Python pyecharts輸出保存圖片代碼實例,文中通過示例代碼介紹的非常詳細,對大家的學習或者工作具有一定的參考學習價值,需要的朋友可以參考下
    2020-11-11

最新評論

精品国内自产拍在线观看