概率圖模型（PGM）：貝葉斯網（Bayesian network）初探

1. 從貝葉斯方法（思想）說起 - 我對世界的看法隨世界變化而隨時變化

用一句話概括貝葉斯方法創始人Thomas Bayes的觀點就是：任何時候，我對世界總有一個主觀的先驗判斷，但是這個判斷會隨着世界的真實變化而隨機修正，我對世界永遠保持開放的態度。

1763年，民間科學家Thomas Bayes發表了一篇名為《An essay towards solving a problem in the doctrine of chances》的論文，

這篇論文發表后，在當時並未產生多少影響，但是在20世紀后，這篇論文逐漸被人們所重視。人們逐漸發現，貝葉斯方法既符合人們日常生活的思考方式，也符合人們認識自然的規律，經過不斷的發展，最終占據統計學領域的半壁江山，與經典統計學分庭抗禮。

讓我們暫時回到Thomas Bayes所處的學術時代18世紀，來一起體會下貝葉斯的思想。

當時數理統計的主流思想是“頻率學派”。所謂頻率學派，舉個例子來說：“有一個袋子，里面裝着若干個白球和黑球（例如3黑2白），請問從袋子中取得白球的概率θ是多少？

頻率派把需要推斷的參數θ看做是固定的未知常數，即概率雖然是未知的，但最起碼是確定的一個值。同時，樣本X是隨機的，所以頻率派重點研究樣本空間，大部分的概率計算都是針對樣本X的分布。

這種觀點看起來確實沒有什么問題，也很合理。但是我們仔細深度思考一下就會發現一個重大問題。

頻率學派之所以能夠得到這種確定性，是因為研究的對象是“簡單可數事件”，例如裝有固定數量球的袋子、只有正反兩面的硬幣、只有6面的標准篩子。但是當研究的問題變得復雜之后，頻率理論就無法處理了。

例如一個朋友創業，你現在需要估計他創業成功的幾率有多大？

這個時候你就無法逐一枚舉出導致他成功或者失敗的所以子原因了（做法有方法？思路清晰？有毅力？能團結周圍的人？和其他競爭對手相比，好多少？....），這是一個連續且不可數的事件空間。

“貝葉斯學派”的觀點和頻率學派則截然相反，貝葉斯學派認為參數是隨機變量，是未知的。而樣本X是固定的，由於樣本是固定的，所以他們重點研究的是參數的分布。

用貝葉斯學派的理論來回答上面創業者評估問題，假如你對這個創業者為人比較了解，你會不由自主的估計他創業成功的幾率可能在80%以上，這是一個先驗的概率估計。隨着公司的運營，你發現公司被運營的非常差，業績也不行，隨即，你對這個創業者的成功估計由80%下調到40%，這是一個后驗概率估計。貝葉斯學派的這種動態的概率估計思維，就是貝葉斯方法。

為了能夠動態地對真實世界進行概率估計，貝葉斯及貝葉斯派提出了一個思考問題的固定模式：

先驗分布 + 樣本信息（觀測結果）  后驗分布

上述思考模式意味着，新觀察到的樣本信息將修正人們以前對事物的認知。換言之，在得到新的樣本信息之前，人們對的認知是先驗分布，在得到新的樣本信息后，人們對的認知為后驗估計。 

筆者思考：

寫到這里的時候，筆者聯想到牛頓宏觀第一定律和愛因斯坦的相對論之間的關系，在宏觀世界，牛頓第一定律是有效的，但是進入高速微觀世界，就只能用更抽象的相對論來概括了。牛頓第一定律只是相對論在宏觀低速世界的一個特例。

同樣，對簡單問題來說，因為事件空間有限可數，所以頻率理論是成立的，但是對於真實世界的復雜問題，事件空間是連續不可數的，就需要用到貝葉斯理論來概括描述了。

Relevant Link:     

《概率圖模型:原理與技術》DaphneKoller 
https://github.com/memect/hao/blob/master/awesome/bayesian-network-python.md

 

2. 一些背景知識

在本章中，我們回顧一些重要的背景材料，這些材料源自概率論、信息論和圖論中的一些關鍵知識，它們都是貝葉斯網的重要概念組成部分。

0x1：條件概率公式

條件概率（又稱后驗概率）就是事件A在另外一個事件B已經發生條件下的發生概率。條件概率表示為P(A|B)，讀作“在B條件下A的概率”，

• 聯合概率：表示兩個事件共同發生的概率。A與B的聯合概率表示為或者
• 邊緣概率（又稱先驗概率）：是某個事件發生的概率。邊緣概率是這樣得到的：在聯合概率中，把最終結果中那些不需要的事件通過合並成它們的全概率，而消去它們（對離散隨機變量用求和得全概率，對連續隨機變量用積分得全概率），這稱為邊緣化（marginalization），比如A的邊緣概率表示為P(A)，B的邊緣概率表示為P(B)。 

0x2：全概率公式和貝葉斯概率公式

基於條件概率公式，我們可以繼續推導出貝葉斯公式。關於這個話題，筆者在另一個文章中進行了詳細的總結，這里不再贅述，為了保持文章自包含性，這里摘取部分筆者認為關鍵的部分。

全概率公式

貝葉斯公式

將其和全概率公式進行對比。會發現以下幾點：

• 貝葉斯公式的分子和全概率公式的分解子項是相同的
• 貝葉斯公式的分母等於全概率公式的累加和

筆者認為，貝葉斯公式可以這么理解：貝葉斯公式表達了一個思想，根據某個出現的結果B可以反推所有可能導致該結果的子原因Ai，而具體根據結果B逆推出每個子原因的條件概率比重，則取決於這個子原因和結果B的聯合概率 ，這是從結果推原因的逆概率計算問題。

0x3：隨機變量的鏈式法則（chain rule）

根據條件概率的定義，我們有：

更一般地，如果a1，....，a_k是事件，那么有如下鏈式分解式，

這個等式稱為條件概率的鏈式法則（chain rule）。 換句話說，我們可以將幾個事件組合的概率表示為關於事件之間的遞歸式依賴關系的乘積。值得注意的是，我們可以用事件的任意順序來擴展這個表達式而保持結果不變。

鏈式法則是一個非常有用的定理，它是貝葉斯網中的核心概念，因為貝葉斯網面對和描述的都是復雜事件。一個復雜事件往往是由大量子事件組成的，不同的子事件間既存在獨立關系，也存在關聯依賴關系。鏈式法則提供了一個理解復雜事件的世界觀，即復雜事件是可以分解的。

0x4：隨機變量的獨立性

關於獨立性，有一個很形象直觀的例子。假設一個思想實驗：如果一個社會中存在一種廣泛的思潮，即重男輕女思潮，90%的家庭都認為生男孩比生女孩好，只要生的不是男孩就是繼續生，直到生了一個女孩為止。 以上為背景，經過10年后，請問社會中男女比率會失衡嗎？

這個問題答案可能有些反直覺，答案是：不管經過多少年，社會中男女比例都保持50%的平衡比例。出現這個現象的原因是，每一胎生男生女都是等概率的，不以父母的意志為轉移。

下面我們定義隨機變量的獨立性概念，

假如或者，則稱事件α和事件β，則稱事件α和事件β在P中獨立，記為。

• 從條件熵的角度解釋獨立性，事件α和事件β在P中獨立，意味着事件事件β對事件α的不確定性不提供額外的信息。
• 從條件概率的角度，我們還可以給出獨立性概念的另一種定義，即分布P滿足，當且僅當。證明過程如下，
  • 使用鏈式法則，，由於α和β獨立，所以，因此，，所以有，

• • 從該公式可以看到，獨立性是一個對稱的概念，即意味着。 

我們需要明白的是，在大多數實際的應用中，隨機變量並不是邊緣獨立的，盡管如此，這種方法的一般形式將是我們求解的基礎。

0x5：隨機變量條件獨立性

關於條件獨立性，筆者通過一個具體的例子來直觀說明。我們假定學生的GPA成績為隨機變量G，被Stanford錄入為隨機變量S，被MIT錄取為隨機變量M。很顯然，我們會得到下列幾個分布判斷：

• 零先驗下條件依賴：”被MIT錄取M“和”被Stanford錄取S“之間並不相互獨立，如果已知一個學生已經被Stanford錄取，那么我們對她被MIT錄取的概率的估計就會提高，因為”被Stanford錄取S“潛在說明了她是一個智商高有前途的學生。換句話說，”被Stanford錄取S“這個事實給”被MIT錄取M“輸入了額外的有效信息熵。
• 給定額外事件下條件獨立：現在假定兩所大學僅根據學生的GPA作決定，並且已知學生的GPA是A。在這種情況下，”被Stanford錄取S“和”被Stanford錄取S“之間就完全獨立了。換句話說，在已知GPA的情況下，”被Stanford錄取S“對”被Stanford錄取S“的推斷沒有提供任何有效的信息熵，這個現象在貝葉斯網中被稱作局部獨立性。用條件概率公式表述上面這句話就是，對於這種情況，我們說在給定Grade條件下，MIT與Stanford條件獨立。

下面我們形式化地定義事件之間條件獨立性的概念，

假如，或者，則稱事件α在給定事件Υ時，在分布P中條件獨立於事件β，記作。

進一步地有，P滿足當且僅當。

現在我們將概念引申到隨機變量的范疇內，關注隨機變量間的條件獨立性。

令X，Y，Z表示隨機變量，如果P滿足，則稱集合X與Y在分布P中條件獨立，集合Z中的變量稱為觀測（observed）變量。觀測（observation）這個詞很能引起我們的思考，我們看到，在不同的觀測下，隨機變量的條件獨立性發生了改變，這正是體現了貝葉斯思想的精髓。

特別的，如果集合Z是空集，可以把記作，並且稱X與Y是邊緣獨立的（marginally independent）。

0x6：查詢一個概率分布

我們基於貝葉斯網進行的一個常見應用是，利用多維隨機變量的聯合概率分布來完成推理查詢。

1. 概率查詢

概率查詢是最常見的查詢，查詢由兩部分組成，

• 證據（條件）：模型中隨機變量的子集E，以及這些變量的實例e
• 查詢變量：網絡中隨機變量的子集Y

概率查詢的任務是計算下式：

即，Y的值y上的后驗概率分布取決於E=e，這個表達式也可以看作是以e為條件作用獲得的分布中Y上的邊緣。

2. 最大后驗概率（MAP）查詢

第二類重要任務是對變量的一些子集找到一個高概率的聯合賦值。這類任務的最簡單變形是MAP查詢（最可能解釋MPE），其目的是找到對所有（非證據）變量最可能的賦值。

如果令，那么在給定證據的條件下，我們的任務是為W中的變量找到最可能的賦值：

其中的表示使得最大的值。

0x7：圖

把圖作為一種數據結構來表示概率分布在貝葉斯網的一個重要概念，我們這節討論一些圖的基本概念。

1. 節點與邊

圖是一個包含節點（node）集與邊（edge）的數據結構。 

假定節點集為，節點對 X_i 與 X_j 由一條有向邊（directed edge）X_i -> X_j 或者一條無向邊（undirected edge）X_i - X_j 連接。因此，邊集是成對節點的集合，圖通過二元關系來存所有邊的連接關系。

因此，圖按照邊集的類型，可以分為：

• 有向圖（directed graph）：通常記作。
• 無向圖（undirected graph）：通常記作，通過忽略邊的方向，有向圖可以轉換為無向圖。

我們用來表示 X_i 與 X_j 經由某種邊連接的情形，這條邊或者是有向的，或者是無向的。

下圖是部分有向圖的一個例子，

其中，

• A是C的唯一父節點
• F，I 是C的子節點
• C的唯一鄰節點是D
• A，D，F，I 均與C相鄰

2. 子圖

很多時候，我們希望只考慮與節點的特定子集相關的部分圖。

1）導出子圖

令，且令，導出子圖（induced subgraph）定義為圖，其中，表示使得的所有邊。例如下圖表示的是導出子圖， 

 

2）完全子圖

如果X中的任意兩個節點均由一條邊連接，那么X上的子圖是完全子圖（complete subgraph）。集合X通常稱為團（clique）。對於節點的任意超集，如果Y不是團，那么X稱為最大團（maximal clique）。

盡管理論上說，節點子集X可以是任意的，但我們一般只對能夠保持圖結構的特定性質的節點集感興趣。

3. 路徑與跡

利用邊的基本符號，可以在圖中對不同類型的長范圍連接進行定義。

1）路徑（path）

對每一個 i=1，...，k-1，如果 X_i -> X_i+1，或者 X_i - X_i+1，那么 X₁，....，X_i 在圖中形成了一條路徑（path）。如果一條路徑中至少存在一個 i，使得 X_i -> X_i+1，則稱該路徑是有向的。

2）跡（trail）

對於每個 i=1，...，k-1，如果，那么 X₁，....，X_i 在圖中形成一條跡（trail）。 

跡是比路徑更弱的約束條件，跡必定是路徑，但路徑不一定成為跡。

例如下圖中，

• A，C，D，E，I 是一條路徑，因此也是一條跡
• A，C，F，G，D 是一條跡，但不是一條路徑

3）連通圖（connected graph）

對於圖中的任意節點 X_i，X_j，如果 X_i 和 X_j 之間存在一條跡，那么該圖是連通圖。

4.  圈（cycle）

1）無圈圖

圖中的一個圈是一條有向路徑 X₁，....，X_k，其中 X₁ = X_k。如果一個圖不包含圈，那么這個圖是無圈圖（acyclic graph）。需要注意的是，在存在有向圈的圖上定義一個一致的概率模型十分困難，因此機器學習中大部分的圖算法都是假設圖是無圈的。

有向無圈圖（directed acyclic graph，DAG）是構成貝葉斯網的基本圖表示。

下圖是一個無圈圖，

但是只要為圖中添加一條無向邊 A-E，就可以獲得一條從A到自身的路徑A，C，D，E，A。

注意，圈是針對路徑的概念而言的，無圈並不代表從節點到自身不會有跡，例如上圖中，包含了 C，D，E，I，C 和 C，D，G，F，C 這兩條跡。

2）部分有向無圈圖（partially directed acyclic graph，PDAG）

一個包含有向邊和無向邊的無圈圖稱為部分有向無圈圖，PDAG上無圈的要求意味着圖可以被分解為一些鏈分支（chain component）的有向圖，其中，每條鏈分支內的節點相互之間僅通過無向邊鏈接。

PDAG的無圈性保證我們可以對分支進行排序，進而使得所有的邊由數值低的分支指向數值高的分支。

3）鏈分支（chain component）

令是上的PDAG，K₁，....，K_i 是的滿足如下條件的不相交划分：

• K_i 的導出子圖不包含有向邊。
• 對 i<j，以及任意一對節點和，X與Y之間的邊僅僅是有向邊 X->Y。 

滿足以上條件時，每一個分支 K_i 稱為一個鏈分支。同時由於其鏈結構，PDAG也稱為鏈圖（chain graph）。

還是對上面的圖例子，我們來看一下鏈分支分解結果，

上圖中共有六個鏈分支：

• {A}
• {B}
• {C，D，E}
• {F，G}
• {H}
• {I}

特別的，當PDAG是一個無向圖時，整個圖形成了一個單鏈分支。相反地，如果PDAG是一個有向圖（無圈的），那么圖中的每一個節點就是其本身的一個鏈分支。

筆者提醒：

鏈圖這個概念非常重要，它和隨機變量獨立性原理一起，構成了貝葉斯網因子分解和局部獨立性原理，也是貝葉斯網能夠從原始的鏈式分解公式得到緊湊結構的理論基礎。

5. 環（loop）

環是一個比圈更弱約束的結構，中的一個環是一條跡 X₁，....，X_k，其中 X₁ = X_k。如果一個圖不包含任何環，那么這個圖是單連通（singly connected）。如果一個節點在單連通圖中只有一個相鄰的節點，那么該節點是一個葉（leaf）節點。

1）森林

如果有向圖的每個節點至多只有一個父節點，那么這個有向圖是森林。

2）樹

如果有向圖的每個節點至多只有一個父節點，且這個有向森林是連通的，那么它是一棵樹。

3）多重樹

一個單連通的有向圖也稱為一棵多重樹（polytree）。注意，多重樹與樹非常不同，例如下圖是多重樹但不是樹，因為有一些節點有不止一個父節點。

多重樹例子

6. 弦（chord）

弦的定義在基於圖結構的推理代價估計中扮演者非常重要的角色。

令 X₁ - X₂ - .... - X_k - X₁ 是圖中的一個環。環中的一條弦（chord）是連接兩個不連貫節點 X_i 與 X_j 的一條邊。對 k>=4，如果任意環 X₁ - X₂ - .... - X_k - X₁ 有一條弦，那么無向圖稱為弦圖（chordal graph）。

值得注意的是，在一個通弦的圖中，最長的”最小環（一個沒有捷徑的環）“是一個三角形。因此，弦圖通常也稱為三角剖分圖（triangulated graph）。

進一步的，如果圖的基礎無向圖是通弦的，則稱圖是通弦的。

Relevant Link:    

《概率圖模型:原理與技術》DaphneKoller 

 

3. 貝葉斯網原理

0x1：什么是貝葉斯網絡

• 貝葉斯網絡(Bayesian network)，又稱信念網絡(Belief Network)，或有向無環圖模型(directed acyclic graphical model)，是一種概率圖模型，於1985年由Judea Pearl首先提出。是一種幫助人們將概率統計應用於復雜領域、進行不確定性推理和數值分析的工具。
• 貝葉斯網是一種系統地描述隨機變量之間關系的語言，它模擬了人類推理過程中因果關系的不確定性。
• 貝葉斯網是概率論與圖論結合的產物，一方面用圖論的語言直觀揭示問題的結構，另一方面按概率論的原則對問題結構加以利用，本質上其網絡拓朴結構是一個有向無環圖(DAG)。
• 許多經典多元概率模型都是貝葉斯網的特例，例如朴素貝葉斯模型、馬爾柯夫鏈、隱馬爾科夫模型、卡爾曼濾波器、條件隨機場等

0x2：為什么需要貝葉斯網絡

在開始介紹貝葉斯網概念之前，我們先來一個討論一下學者們提出貝葉斯網的學術動機，筆者將其總結為下面幾句話：

• 不確定性推理
• 高維稀疏表征問題：利用條件獨立降低復雜度
• 局部概率模型（local probability model）：貝葉斯網的因子分解
• 有限狀態轉移的動態貝葉斯過程：因果推理與證據推理

1. 不確定性推理

所謂不確定性推理，是指基於概率論對真實世界進行抽象建模以及數值分析，使用概率分布而不是精確數值的方式進行有效預測。常規的不確定性推理過程如下，

• 把問題用一組隨機變量 X = {X1, X2, . . . , Xn} 來刻畫
• 把關於問題的知識表示為一個聯合概率分布 P(X)
• 按概率論原則進行推理計算

2. 利用條件獨立降低復雜度

根據全概率公式和鏈式法則我們知道，每個隨機事件都可以分解為大量的子事件原因的累乘，且子事件原因之間還存在大量的遞歸依賴，

這樣，假設我們的目標是在某個隨機變量的集合上表示聯合分布P，即使在變量為二值（binary-valued）這樣最簡單的情況下，n個變量的聯合分布也需要個不同賦值的概率。當n較大時，這種聯合分布的顯示表示和樣本學習都會變得非常難處理。

• 從計算的角度，操作的成本非常昂貴，並且也很難在內存中存儲
• 從機器認知的角度，我們不可能從專家那里獲得如此多的數值。而且，這些數值非常小並且不能與人們可以合理思考的事件一一對應
• 從統計的角度，如果希望從數據中學習出分布，需要有規模非常龐大的數據才能可靠地估計出這么多的自由參數，否則會遇到欠擬合問題

但是實際上，正如我們前面介紹舉例的，對於一個好學生來說，考上Stanford和考上MIT這兩件事之間是獨立的，條件概率和邊緣概率是相同的。因此，貝葉斯網絡基於隨機變量間獨立性原理，對去全概率公式進行了有效的化簡，得到了一個更緊湊的高維分布。這對我們建模復雜性事件提供了條件，例如下圖，

學生示例的貝葉斯網圖

可以以兩種不同的方式來看待這個圖：

• 它是提供了以因子分解的方式緊湊表示聯合分布骨架的數據結構
• 它是關於一個分布的一系列條件獨立性假設的緊湊表示

在嚴格意義上，這兩種觀點是等價的。

注意到，貝葉斯網用邊表示直接依賴關系，如果在真實世界中A只依賴於B，則只應該存在一條從B通向A的有向邊。節點只直接依賴其父節點，是貝葉斯網語義的核心。

3. 基於局部概率模型表示成對的因果推理

為了說明這個概念，我們繼續引用學生示例，

這個問題包含五個隨機變量：

• 學生的智商（I）
• 課程的難度（D）
• 考試的成績（G）
• 學生的SAT成績（S）
• 推薦信的質量（L）

我們知道，貝葉斯網可以用有向圖來表示，其節點表示隨機變量，而邊則表示一個變量對另一個變量的直接影響（現實世界運行方法的直觀假設）。我們可以將這個圖看成是對一個由自然執行的生成式采樣過程進行編碼，其中，每個變量的值指取決於其父節點的分布自然選擇。換句話說，每個變量是其父節點的一個隨機函數（局部條件概率性質）。

換句話來說，局部條件概率性質來源於問題域本身的： 

• 學生的成績不僅取決於他的智商（I），同時也取決於課程的難度（D）
• 學生請求教授為其學推薦信，不過教授健忘並且從來記不住學生的名字，所以教授通過查看學生的成績，並依據這僅有的信息來為其寫推薦信，推薦信的好壞用L表示，隨機地依賴於考試的成績
• SAT成績僅僅取決於學生的智商，在智商（I）已知的情況下，SAT成績（S）和課程成績（G）是獨立的

我們將貝葉斯網中包含的這種局部的條件概率結構稱之為局部概率模型（local probability model），用於表示每個變量對其父節點依賴的本質。 

對於這個問題，條件概率分布的一個可能選擇如下圖所示，這個網絡結構與條件概率分布一起組成了貝葉斯網，

帶有條件概率分布的學生示例的貝葉斯網

具體的，考慮空間中的一些特定狀態，比如，直觀上，這個事件的概率可以通過計算構成它的基本事件的條件概率乘積獲得，即：

一般地，對於聯合概率空間中的任何狀態，可以采用相同的過程，即：

上述本質上是一個概率查詢（推理）過程。 

4. 有限狀態轉移的動態貝葉斯過程

1）因果推理

所謂的“推理”，顧名思義是一個根據當前獲得的證據動態變化條件概率的過程。現在考慮貝葉斯網，並且觀察當證據獲得時各種事情的概率如何動態地發生變化。

以某個學生George為例，用獲得的模型對George進行推理。我們可能想要詢問George有多大可能從教授Econ101的教授那里獲得一封好的推薦信（l¹）。現在我們開始動態貝葉斯推理：

• 在對George和Econ101的其他情況一無所知的條件下，P(l¹)=50.2%。更精確地，令為根據前述貝葉斯網（BN）定義的聯合分布，那么。
• 現在我們了解到George不太聰明（i⁰），現在George從教授Econ101課程的教授那里獲得一封好的推薦信的概率降低到大約38.9%，即。

• 現在我們進一步得知Econ101是一門簡單的課程（d⁰），此時，George從教授那里獲得一封好的推薦信的概率是。

• 注意到上面兩步的動態推理，都是基於我們不知道George的成績Grade前提下，即head-to-head的關鍵節點未知，如果Grade已知，智商和課程難度就會獨立，互不影響。

這類可以預測各種因素（如George的智商）“順流而下”影響的概率查詢，我們稱之為因果推理（causal reasoning）或者預測（prediction）。

2）證據推理

現在考慮Acme咨詢公司的一個招聘官，他試圖貝葉斯網決定是否雇佣George，我們開始動態貝葉斯推理：

• 先驗知識是，招聘官相信George有30%的可能有高智商。
• 之后，他拿到了George在Econ101上的成績記錄為C（g³），此時，招聘官對George有高智商的概率降至7.9%，即，同時，課程困難度的概率也從40%上升到62.9%，因此我們有理由相信是因為課程太難導致George沒有取得高分。

• 現在，假定招聘官不小心丟失了George的成績單，關於George的信息只剩下來自其Econ101課程教授的推薦信，不過這封信不是一個好的推薦信，此時，George有高智商的概率依然會降低，但只降至14%，即。我們看到，節點距離越近，信息傳遞的影響就越大，這也符合馬爾科夫過程的熵傳遞原理。

•  同時我們注意到，如果招聘官既有成績單又有推薦信，George有高智商的概率與招聘官只有成績單時的概率是一樣的，。原因是在head-to-tail結構下，如果關鍵節點已知，則head與tail之間是條件獨立的。

• 最后，George向招聘官提交了他的SAT成績，他的SAT成績很高，George有高智商的概率戲劇性地從7.9%又上升到了57.8%，，即。產生這個現象的原因是，與不好的考試成績相比，SAT高分更有價值，因為低智商的學生在SAT考試中獲得高分的可能性非常低，遠遠低於Econ101課程考試，也可以通俗地說“SAT考試成績更硬”。同時，因為高智商的學生在困難的課程中更可能獲得C，所有現在Econ101是一門困難的課程的概率已經從62.9%提高了76%左右。 

這類從結果推理原因的“逆流而上”的概率查詢，我們稱之為證據推理（evidential reasoning）或解釋（explanation）。 

3）因果間推理一般模式

我們知道，課程難度D，課程成績G，智商I，它們之間是head-to-head的拓朴結構。

如果George在Econ101中得到的成績是B，那么。另一方面，如果Econ101比較難，那么。

實際上我們已經通過課程的難度解釋消除（explaining away）了成績低這一事件帶來的影響。解釋消除是稱為因果間推理（intercausal reasoning）的一般推理模式的一個示例。其中相同結果的不同原因之間可以相互影響。 在這里，因為Grade已知，所以D和I之間的信息影響流是通的，影響力可以沿着跡傳遞。

這類推理在人類日常的推理中是非常普遍的，例如，當發燒、咽喉痛，並且擔心可能有單核細胞增多症時，如果告知只是得了流感，那么得單核細胞增多症的概率就會被消解下降。因為，患有流感為出現的針狀提供了另一種解釋，從而很大程度上降低了患有單核細胞增多症的概率。

筆者思考：

對上述這個理論可以有一個更加形象的理解方式，即河流類比，我們將新的證據帶來的信息熵比作河流里的水，學生在某門課程上取得了不好的成績就是一個有效信息熵，這個信息順着G-I之間流動，使得高智商的后驗概率下降。現在加入一個新的信息源是Econ101課程比較難，這個新的信息源就相當於一個新的水源，沿着D-G-I這個河道傳播，和原來的G-I一起，共同綜合決定了I的后驗概率。

同樣的，如果G是未知的，那就相當於這個河流中間是斷開的，無論D如何都無法影響到I，I的后驗概率只取決於G。

4）有限狀態轉移

另外一方面，貝葉斯網必須是有向無圈圖。其中節點代表隨機變量，節點間的邊代表變量之間的直接依賴關系。

本質上因為貝葉斯網是對真實世界復雜性事件的抽象建模，所以不允許出現圈（即循環遞歸依賴）。簡單來說就是不允許出現下列這種場景，

• A事件的發生依賴於B
• B事件的發生依賴於C
• C事件的發生依賴於D
• D事件的發生依賴於A

如果出現循環依賴，因果分解鏈式公式就會無限膨脹，也就無法用計算機處理。

0x3：貝葉斯網絡的定義

貝葉斯網絡的有向無環圖中的節點表示隨機變量，它們可以是可觀察到的變量，或隱變量、未知參數等。

有因果關系（或非條件獨立）的變量或命題，用箭頭來連接。若兩個節點間以一個單箭頭連接在一起，表示其中一個節點是“因(parents)”，另一個是“果(children)”，兩節點就會產生一個條件概率值。

例如，假設節點E直接影響到節點H，即E->H，則用從E指向H的箭頭建立結點E到結點H的有向弧(E,H)，權值(即連接強度)用條件概率P(H|E)來表示，如下圖所示：

一般地，貝葉斯網結構是其節點代表隨機變量X₁，....，X_n的一個有向無圈圖（DAG）。令表示X_i在中的父節點，表示X_i在圖中的非后代節點變量。因此表示了如下稱為局部獨立性（local independencies）的條件獨立性假設，並且記為。

對每個變量X_i，都有：

換句話說，局部獨立性表明，在給定父節點條件下，每個節點X_i與其非后代節點條件獨立。貝葉斯網圖的形式化語義是一系列獨立性斷言。 

如下圖所示，便是一個簡單的貝葉斯網絡：

因為a導致b，a和b導致c，所以有下面累乘公式，

 

0x4：貝葉斯網絡的典型結構形式

3節點對應的三角結構是貝葉斯網信息流動的最基本結構，這節我們來介紹幾種經典的信息流動結構。

1. head-to-head（共同的作用）

根據上圖我們有：P(a,b,c) = P(a) * P(b) * P(c|a,b)成立，化簡后可得：

這種結構代表的意思是：

• c未知時：a、b被阻斷(blocked)，是獨立的，稱之為head-to-head條件獨立。
• c已知時（或已知c的后代）：a和b不獨立，a的信息可以沿着c流向b。

我們可以形象地理解為，將abc想象成河流，兩條小河（a、b）流入一條大河（c），當c是未知時，a->c<-b之間的流動是阻斷的，那么影響將無法沿着“a->c<-b”流動。形如“a->c<-b”的這種結構也稱為一個v-結構（v-structure）。

2. tail-to-tail（共同的原因）

這種結構代表的意思是：

• c未知時：有：P(a,b,c)=P(c) * P(a|c) * P(b|c)，此時，沒法得出P(a,b) = P(a)P(b)，即c未知時，a、b不獨立，a可以通過c來影響b。
• c已知時：有：P(a,b|c)=P(a,b,c) / P(c)，然后將P(a,b,c)=P(c) * P(a|c) * P(b|c)帶入式子中，得到：P(a,b|c)=P(a,b,c)/P(c) = P(c)*P(a|c)*P(b|c) / P(c) = P(a|c)*P(b|c)，即c已知時，a，b被阻斷(blocked)，a、b是獨立的，稱之為tail-to-tail條件獨立。

我們可以形象地理解為，一條大河到某處分叉成兩條支流，在c給定的條件下，a，b是兩條獨立的支流，彼此是獨立的。

3. head-to-tail（因果跡）

這種結構代表的意思是：

• c未知時：有：P(a,b,c)=P(a)*P(c|a)*P(b|c)，但無法推出P(a,b) = P(a)P(b)，即c未知時，a、b不獨立。
• c已知時：有：P(a,b|c)=P(a,b,c)/P(c)，且根據P(a,c) = P(a)*P(c|a) = P(c)*P(a|c)，可化簡得到：P(a,b|c) = P(a,b,c)/P(c) = P(a)*P(c|a)*P(b|c)/P(c) = P(a,c) * P(b|c) / P(c) = P(a|c) * P(b|c)，即在c給定的條件下，a，b被阻斷(blocked)，是獨立的，稱之為head-to-tail條件獨立。 

4. tail-to-head（證據跡）

這種結構代表的意思是： 

• c未知時：影響可以經過c從a流向b，a可以作為證據影響b的后驗概率
• c已知時：在c給定的條件下，a，b被阻斷，是獨立的。

0x5：條件概率鏈式法則與局部獨立性的統一

我們知道，貝葉斯網的形式化語義是一系列獨立性斷言，另一方面，貝葉斯網也是由一系列條件概率分布所組成的圖，這個圖通過鏈式法則為貝葉斯網定義了一個聯合分布。

在本節中，我們將討論它們二者的等價性。分布P滿足於圖相關的局部獨立性，當且僅當P可以由與圖相關的一系列條件概率分布表示。

我們從形式化基本概念開始。

1. I-Map

首先我們定義與分布P相關的獨立性的集合。

令P是上的分布，定義為在P中成立的形如的獨立性斷言的集合。

現在可以將陳述“P滿足與相關的局部獨立性”簡寫為。在這種情況下，可以稱是P的一個I-map(獨立圖)。

更一般的定義如下：

令是與獨立性集合相關聯的任意一個圖對象，如果對獨立性集合，，則稱是一個I-map。 

2. I-map到因子分解

一個貝葉斯網蘊含了一系列條件獨立性假設，這是我們能夠用緊湊因子分解表示貝葉斯網的關鍵。一個特殊的例子就是朴素貝葉斯分解公式。

考慮學生貝葉斯網的例子，

先考慮任意分布P，我們將對聯合分布進行分解並且將其分解成一些局部概率模型：

 

這種轉換不依賴任何假設，它對任何聯合分布P都適用。然而，因為等式右邊的因子分解中的條件概率沒有充分利用我們已知的獨立依賴性信息，因此非常不緊湊。

現在，我們應用由貝葉斯網導出的條件獨立性假設，假定關於分布P是一個I-map，我們得到下式：

• 
• 
• 因為，所以有

綜上，聯合分布可化簡為：

從上式可以看出，聯合分布中的任何表值都可以通過每個變量所對應的因子的乘積來計算。每個因子表示表示一個變量在網絡中給定父節點時的條件概率。在是I-map的條件下，這種因子分解可以應用於任何分布P。

用形式化的語言來表達就是：

令未定義在變量X₁，...，X_n上的一個貝葉斯網，假如P可以表示為如下乘積：

則稱分布P是關於圖的在同一空間上的因子分解。這個等式稱為貝葉斯網的鏈式法則。單個因子稱為條件概率分布(CPD)或局部概率模型。

貝葉斯網結構所蘊含的條件獨立性假設允許我們將是一個I-map的分布P分解成一系列較小的條件概率分布。

0x6：d-分離概念

上一節我們討論了聯合分布的因子分解和局部獨立性的統計概念。這個小節我們將其推廣到整個貝葉斯網中，首先我們定義有效跡的概念。

1. 有效跡

考慮跡比較長的一般情況。直觀上，對於從X₁流動到X_n的影響，它需要流經這條跡上的單個節點。換句話說，如果每條雙邊跡都允許影響流過，那么X₁可以影響到X_n。

對這一直觀認識形式化定義如下：

令是一個貝葉斯網結構，且是中的一條跡。令Z是觀測變量的一個子集，在給定Z的條件下，如果滿足下列條件：

• 一旦有一個v-結構，則X_i或其一個后代在Z中
• 跡上的其他節點都不在Z中

那么稱跡是有效跡。

2. d-分離（d-separation）

一般的，如果節點之間存在任何允許影響流過的跡，那么一個節點可以影響另一個節點。基於這個直觀認識，可以定義d-分離的概念，這個概念為我們提供了在有向圖的節點之間實現分離的能力。

形式化定義如下：

令X，Y，Z是圖的三個節點集，在給定Z的的條件下，加入任意節點與之間不存在有效跡，那么X與Y在給定Z時是d-分離的，記作。與d-分離相對應的獨立性的集合用表示：

這個集合也稱為全局馬爾科夫獨立性（global Markov independencies）集。中的獨立性恰好是那些可以保證在的每個分布中都成立的獨立性。 

0x7：構建貝葉斯網結構的一般化方法

從上一章的條件概率鏈式公式化簡過程我們可以看到，對於同一個獨立性集合，與其相符的結構有很多。一個好的方法是選擇能夠反映真實世界的因果序列和依賴關系的結構，目的是使原因成為結果的父節點。

例如，在對汽車保險相關的貝葉斯網構建中，以前的事故（A）往往是好司機（G）的父節點，即“A->G”，這與保險公司對汽車保險這一問題的認知一致。同時，由於司機不好（B）才造成了以前（以及將來的）事故，即“B->A”。

 

一般來說，建模時我們可能會選擇很多比較弱的影響，但如果將它們全部考慮進去，網絡會變得非常復雜。這樣的網絡從表示的角度看是有問題的，它們難於理解、難於修正和學習。此外，由於貝葉斯網中的推理非常強地依賴於它們的連接功能，因此，增加這樣的邊無疑會使網絡訓練和預測的代價變得很高。

0x8：貝葉斯網學習（結構構件和參數估計）

貝葉斯網的學習過程分為兩部分：

• 結構學習：發現變量之間的圖關系，構建貝葉斯網結構的一個一般方法是后向構建過程，該過程從我們感興趣的一個變量開始，比如肺癌。然后嘗試找到這個變量的先驗因素，並將其添加為這個變量的父節點，這個過程稱為對話擴展。例如，我們可以用這個過程來得出如下結論：
  • 肺癌應該把吸煙作為其父節點
  • 吸煙（可能，但不明顯）應該把性別和年齡作為其父節點
• 參數學習：決定變量之間互相關聯的量化關系，即估計貝葉斯網絡的 CPD（條件概率分布） 表格中的數值。操作方法很簡單，就是計算訓練數據中事件發生的次數。例如要估計 p(SAT=s1|Intelligence=i1)，我們只需要計算 SAT=s1 且 Intelligence = i1 的數據點在 Intelligence = i1 的數據點總量中所占的比例。

從數據中學習貝葉斯網絡 

0x9：貝葉斯網推理

所謂推理，就是要根據我們已知的信息做出預測。我們可以使用推理來解答一些問題：

1. 邊際推理（marginal inference）

尋找一個特定變量的概率分布。比如，給定一個帶有變量 A、B、C 和 D 的圖，其中 A 取值 1、2 和 3，求 p(A=1)、p(A=2) 和 p(A=3)。

2. 后驗推理（posterior inference）

給定某些顯變量 v_E（E 表示證據（evidence）），其取值為 e，求某些隱藏變量 v_H 的后驗分布 p(v_H|v_E=e)。

3. 最大后驗（MAP）推理（maximum-a-posteriori inference）

給定某些顯變量 v_E，其取值為 e，求使其它變量 v_H 有最高概率的配置。

Relevant Link:    

《概率圖模型:原理與技術》DaphneKoller 
http://staff.ustc.edu.cn/~jianmin/lecture/ai/bnet_slides.pdf
https://zhuanlan.zhihu.com/p/30139208
https://www.zhihu.com/question/28006799
https://blog.csdn.net/v_JULY_v/article/details/40984699
http://staff.ustc.edu.cn/~jianmin/lecture/ai/bnet_slides.pdf

 

4. 貝葉斯網的一個特殊特例 - 馬爾柯夫鏈

關於馬爾柯夫鏈，筆者在另一篇文章有詳細的討論。簡單來說，馬爾柯夫鏈是一種基於“馬爾科夫假設（marko assumption）”的動態貝葉斯網。

從網絡結構的角度來說，馬爾柯夫鏈是一個head-to-tail鏈式網絡，如下圖所示：

我們已經知道，在xi給定的條件下，xi+1的分布和x1,x2…xi-1條件獨立。這意味着：xi+1的分布狀態只和xi有關，和其他變量條件獨立。

通俗地說，當前狀態只跟上一狀態有關，跟上上或上上之前的狀態無關。這種順次演變的隨機過程，就叫做馬爾科夫鏈（Markov chain）。且有：

可以看到，馬爾柯夫鏈是一個特殊的貝葉斯網結構。

Relevant Link:    

https://www.cnblogs.com/LittleHann/p/7609060.html#_lab2_2_1

 

5. 貝葉斯網的一個特殊特例 - 朴素貝葉斯模型

這章我們來描述一個最簡單的例子，在這里例子中，條件參數化與條件獨立性假設被結合在一起，目的是對高維概率分布產生非常緊湊的表示。

0x1：場景描述

還是前面公司雇佣大學畢業生的例子，現在假定公司除了智商I和SAT成績S之外，同時還知道學生某些課程的成績G，在這種情況下，概率空間是定義在三個相關隨機變量 I、S、G 上的聯合分布。假定 I 與 S 如前所述，並且G可能取三個值g¹，g²，g³，分別代表成績A，B和C，那么聯合分布有12個表值。

但是在這個例子中，我們發現了潛在的冗余結構，分布P滿足條件獨立性執行，即如果我們知道一個學生是高智商的，那么他在SAT考試中獲得的高分並不能為我們帶來任何有關他課程成績的信息，形式化的表達為：

更一般的情況下，可以表達為：

 

注意到，這種獨立性聲明只有當假定學生的智商是其課程與SAT得分的唯一原因時才成立（tail-to-tail結構）。

基於條件獨立性假設，我們可以得到一種更加緊湊的因果推理表達方式，

  

這樣，聯合分布P(I，S，G)可以分解為三個條件概率分布（CPD）的乘積。

作為一個整體，這3個條件概率分布詳細說明了聯合分布在條件獨立性假設下的因果表示，例如：

我們注意到，除了更緊湊的表達方式之外，基於條件獨立性的因果表達方法的另一個好處是，模塊性。當增加一個新的變量G時，聯合分布在整體上會發生變化。而在因子表示中，我們可以重復使用變量I和S的局部概率模型，並且只需增加對新變量的條件概率P(G|I)即可。這個性質對於模擬真實系統非常有價值。

0x2：朴素貝葉斯模型（naive Bayes） 

上個小節討論的例子被稱之為朴素貝葉斯模型。

朴素貝葉斯模型假設所有的實例屬於若干兩兩互斥且包含所有事例情況的類（class）中的一個。因此，存在一個在某個集合{c¹，...，c^k}中取值的類變量C。在這個例子中，類變量是學生的智商I，並且存在且只存在事例的兩個類（高智商和低智商）的學生。

模型還包括一定數量的、可以觀測到其值的特征（feature）X₁，....，X_k。朴素貝葉斯假設（naive Bayes assumption）是在給定事例的類的條件下，這些特征條件獨立。換言之，在事例的每個類內，不同的性質可以獨立地確定，形式化的表達為。

從網絡結構上來說，朴素貝葉斯是一種典型的tail-to-tail的貝葉斯網結構，整個貝葉斯網都是由這種結構組成，

朴素貝葉斯模型的貝葉斯網

基於這些獨立性假設，朴素貝葉斯模型的因子分解（factorization）可以表示如下，

0x3：朴素貝葉斯模型應用之困

盡管具有很強的假設，由於簡單並且只需要少量參數，朴素貝葉斯模型經常在實踐中作為“分類器”得以應用，例如基於朴素貝葉斯模型的垃圾郵件分類器。

這個模型最早被應用於醫療診斷，其中，類變量的不同值用於表示患者可能患的不同疾病。證據變量用於表示不同症狀、化驗結果等。在簡單的疾病診斷上，朴素貝葉斯模型確實發揮了很好的作用，甚至比人類專家的診斷結果都要好。但是在更深度的應用中，醫生發現，對於更復雜（由多種致病原因和症狀共同表現）的疾病，模型表現的並不好。

數據科學家經過分析認為，出現這種現象的原因在於：模型做了集中通常並不真實的強假設，例如：

一個患者至多可能患一種疾病

在已知患者的疾病條件下，不同症狀的出現與否，不同化驗結果，之間是互相獨立的

這種模型可用於醫學診斷是因為少量可解釋的參數易於由專家獲得，早期的機器輔助醫療診斷系統正式建立在這一技術之上。

但是，之后更加深入的實踐表明，構建這種模型的強假設降低了模型診斷的准確性，尤其是“過度計算”某些特定的證據，該模型很容易過高估計某些方面特征的影響。

例如，“高血壓”和“肥胖症”是心臟病的兩個硬指標，但是，這兩個症狀之間相關度很高，高血壓一般就伴隨着肥胖症。在使用朴素貝葉斯公式計算的時候，由於乘法項的緣故，關於這方面的證據因子就會被重復計算，如下式：

P(心臟病 | 高血壓，肥胖症) = P(高血壓 | 心臟病) * P(高血壓 | 肥胖症) / P(高血壓，肥胖症)

由於“高血壓”和“肥胖症”之間存在較強相關性的緣故，我們可以很容易想象，分子乘積增加的比率是大於分母聯合分布增加的比率的。因此，當分子項繼續增加的時候，最終的后驗概率就會不斷增大。但是因為新增的特征項並沒有提供新的信息，后驗概率的這種增大變化反而降低了模型的預測性能。

實際上，在實踐中人們發現，朴素貝葉斯模型的診斷性能會隨着特征的增加而降低，這種降低常常歸因於違背了強條件獨立性假設。

Relevant Link:     

https://www.cnblogs.com/LittleHann/p/7199242.html#_label6

 

6. 基於貝葉斯網解決三門問題（蒙提霍爾問題）

本章來自雲棲社區一篇很棒的文章。

主持人會向挑戰者展示三扇關着的門，其中一扇門之后有一輛車，其它門后則有一些無價值的東西。挑戰者並不知道哪扇門背后是車。

現在挑戰者可以先選擇一扇門。然后，主持人會打開剩下的兩扇門中沒有車的一扇。現在，挑戰者可以選擇是否更換選擇的門，但是問題來了，我們應該更換嗎？

直覺上看，主持人似乎並沒有透露任何信息。事實證明這種直覺並不完全正確。讓我們使用我們的新工具"貝葉斯網"來分析這個問題。

我們首先先定義一些隨機變量：

• D：背后有車的門：可取值為 1、2 或 3；D是未被觀察到的隨機變量。
• F：你的第一個選擇：可取值為 1、2 或 3；D是已被觀察到的隨機變量。
• H：主持人打開的門：可取值為 1、2 或 3；主持人打開其中一扇門之前，H 都是未被觀察到的。
• I：F 是否是 D？：可取值為 0、1；I是未被觀察到的隨機變量。

我們建立初始狀態貝葉斯網結構（即挑戰者做出的選擇，但是主持人還未打開一扇門）：

注意到：

• D 和 F 是相互獨立的，因為此時I和H都是未知的（參考文章前面談到的head-to-head獨立結構）
• I 依賴於 D 和 F
• 主持人選擇的門 H 也取決於 D 和 F
• 目前我們對 D 還一無所知，因為F是和D是獨立的，因此知道F並不能對D有任何信息影響（這與學生網絡的結構類似，即知道學生的智力水平不能讓你獲得有關課程難度的任何信息。） 

現在，主持人選擇了門 H 並打開了它。所以現在 H 已被觀察到。 

現在貝葉斯網的條件概率情況發生了變化：

• I<-D->H是一種tail-to-tail結構，D是未知的，所以 H 的信息影響力可以通過D傳遞到 I，因此已知 H 提供了對 I 的一些信息

也就是說，當支持人打開了剩下的兩扇門中沒有車的一扇后，挑戰者理性的做法是改變自己最初的決策。

下面用CPD表格來形式化地分析這個游戲過程：

• D 和 F 都是先驗概率，背后有車的門可能是這些門中的任何一扇且概率相等，我們選擇其中一扇的概率是一樣的。
• I 依賴於 D 和 F 的先驗：
  • 當 D，F 一樣時，I=1
  • 當 D，F 不一樣時，I=0
• H 依賴於 D 和 F 的先驗：
  • D，F 一樣時，那么主持人從另外兩扇門選擇一扇門的概率一樣
  • 如果 D 和 F 不一樣，那么主持人就選擇第三扇門

現在，讓我們假設我們已經選擇了一扇門。也就是說現在已經觀察到 F，假設 F=1。那么給定 F 時，I 和 D 的條件概率是多少？

基於上面CPD表格計算邊緣概率得：

到目前為止，我們選對了門的概率都是三分之一，汽車仍然有可能在任何一扇門之后且概率相等。

現在，主持人打開了 F 之外的另一扇門，所以我們觀察到了 H。假設 H=2。讓我們再計算給定了 F 和 H 時 I 和 D 的條件概率。

概率計算結果如下：

可以看到，我們第一個選擇正確的概率仍然是三分之一，我們的直覺也是如此。但是，現在車在第 3 扇門后的概率不再是三分之一，而是三分之二了。

所以如果我們更換選擇，那么我們得到車的概率是三分之二；如果我們不換，我們得到車的概率是三分之一。 

示例代碼如下：

import math
from pomegranate import *

# We'll create the discrete distribution for our friend first.
friend = DiscreteDistribution( { True: 0.5, False: 0.5 } )

# The emissions for our guest are completely random.
guest = ConditionalProbabilityTable(
    [
        [True, 'A', 0.50],
        [True, 'B', 0.25],
        [True, 'C', 0.25],
        [False, 'A', 0.0],
        [False, 'B', 0.7],
        [False, 'C', 0.3]
    ],
    [friend]
)

# Then the distribution for the remaining cars.
remaining = DiscreteDistribution({0: 0.1, 1: 0.7, 2: 0.2, })

# The probability of whether the prize is randomized is dependent on the number of remaining cars.
randomize = ConditionalProbabilityTable(
    [
        [0, True, 0.05],
        [0, False, 0.95],
        [1, True, 0.8],
        [1, False, 0.2],
        [2, True, 0.5],
        [2, False, 0.5]
    ],
    [remaining]
)

# Now the conditional probability table for the prize. This is dependent on the guest's friend and whether or not it is randomized.
prize = ConditionalProbabilityTable(
    [
        [True, True, 'A', 0.3],
        [True, True, 'B', 0.4],
        [True, True, 'C', 0.3],
        [True, False, 'A', 0.2],
        [True, False, 'B', 0.4],
        [True, False, 'C', 0.4],
        [False, True, 'A', 0.1],
        [False, True, 'B', 0.9],
        [False, True, 'C', 0.0],
        [False, False, 'A', 0.0],
        [False, False, 'B', 0.4],
        [False, False, 'C', 0.6]
    ],
    [randomize, friend]
)

# Finally we can create the conditional probability table for our Monty. This is dependent on the guest and the prize.
monty = ConditionalProbabilityTable(
    [
        ['A', 'A', 'A', 0.0],
        ['A', 'A', 'B', 0.5],
        ['A', 'A', 'C', 0.5],
        ['A', 'B', 'A', 0.0],
        ['A', 'B', 'B', 0.0],
        ['A', 'B', 'C', 1.0],
        ['A', 'C', 'A', 0.0],
        ['A', 'C', 'B', 1.0],
        ['A', 'C', 'C', 0.0],
        ['B', 'A', 'A', 0.0],
        ['B', 'A', 'B', 0.0],
        ['B', 'A', 'C', 1.0],
        ['B', 'B', 'A', 0.5],
        ['B', 'B', 'B', 0.0],
        ['B', 'B', 'C', 0.5],
        ['B', 'C', 'A', 1.0],
        ['B', 'C', 'B', 0.0],
        ['B', 'C', 'C', 0.0],
        ['C', 'A', 'A', 0.0],
        ['C', 'A', 'B', 1.0],
        ['C', 'A', 'C', 0.0],
        ['C', 'B', 'A', 1.0],
        ['C', 'B', 'B', 0.0],
        ['C', 'B', 'C', 0.0],
        ['C', 'C', 'A', 0.5],
        ['C', 'C', 'B', 0.5],
        ['C', 'C', 'C', 0.0]
    ],
    [guest, prize]
)

# Now we can create the states for our bayesian network.
s0 = State(friend, name="friend")
s1 = State(guest, name="guest")
s2 = State(prize, name="prize")
s3 = State(monty, name="monty")
s4 = State(remaining, name="remaining")
s5 = State(randomize, name="randomize")

# Now we'll create our bayesian network with an instance of BayesianNetwork, then add the possible states.
network = BayesianNetwork("test")
network.add_states(s0, s1, s2, s3, s4, s5)

# Then the possible transitions.
network.add_transition(s0, s1)
network.add_transition(s1, s3)
network.add_transition(s2, s3)
network.add_transition(s4, s5)
network.add_transition(s5, s2)
network.add_transition(s0, s2)

# With a "bake" to finalize the structure of our network.
network.bake()

# Now let's create our network from the following data.
data = [
    [True,  'A', 'A', 'C', 1, True],
    [True,  'A', 'A', 'C', 0, True],
    [False, 'A', 'A', 'B', 1, False],
    [False, 'A', 'A', 'A', 2, False],
    [False, 'A', 'A', 'C', 1, False],
    [False, 'B', 'B', 'B', 2, False],
    [False, 'B', 'B', 'C', 0, False],
    [True,  'C', 'C', 'A', 2, True],
    [True,  'C', 'C', 'C', 1, False],
    [True,  'C', 'C', 'C', 0, False],
    [True,  'C', 'C', 'C', 2, True],
    [True,  'C', 'B', 'A', 1, False]
]

network.fit(data)

# We can see the results below. Lets look at the distribution for our Friend first.
print friend

# Then our Guest.
print guest

# Now the remaining cars.
print remaining

# And the probability the prize is randomized.
print randomize

# Now the distribution of the Prize.
print prize

# And finally our Monty.
print monty

這里介紹一個pomegranate庫，它抽象了概率圖模型的底層細節，我們可以方便地基於API進行上層應用建模。

# -*- coding: utf-8 -*-

from pomegranate import *


if __name__ == '__main__':
    guest = DiscreteDistribution({'A': 1./3, 'B': 1./3, 'C': 1./3})
    prize = DiscreteDistribution({'A': 1./3, 'B': 1./3, 'C': 1./3})
    monty = ConditionalProbabilityTable(
            [['A', 'A', 'A', 0.0],
             ['A', 'A', 'B', 0.5],
             ['A', 'A', 'C', 0.5],
             ['A', 'B', 'A', 0.0],
             ['A', 'B', 'B', 0.0],
             ['A', 'B', 'C', 1.0],
             ['A', 'C', 'A', 0.0],
             ['A', 'C', 'B', 1.0],
             ['A', 'C', 'C', 0.0],
             ['B', 'A', 'A', 0.0],
             ['B', 'A', 'B', 0.0],
             ['B', 'A', 'C', 1.0],
             ['B', 'B', 'A', 0.5],
             ['B', 'B', 'B', 0.0],
             ['B', 'B', 'C', 0.5],
             ['B', 'C', 'A', 1.0],
             ['B', 'C', 'B', 0.0],
             ['B', 'C', 'C', 0.0],
             ['C', 'A', 'A', 0.0],
             ['C', 'A', 'B', 1.0],
             ['C', 'A', 'C', 0.0],
             ['C', 'B', 'A', 1.0],
             ['C', 'B', 'B', 0.0],
             ['C', 'B', 'C', 0.0],
             ['C', 'C', 'A', 0.5],
             ['C', 'C', 'B', 0.5],
             ['C', 'C', 'C', 0.0]], [guest, prize])

    s1 = Node(guest, name="guest")
    s2 = Node(prize, name="prize")
    s3 = Node(monty, name="monty")

    model = BayesianNetwork("Monty Hall Problem")
    model.add_states(s1, s2, s3)
    model.add_edge(s1, s3)
    model.add_edge(s2, s3)
    model.bake()

    print model.probability([['A', 'A', 'A'],
                             ['A', 'A', 'B'],
                             ['C', 'C', 'B']])

    print model.predict([['A', 'B', None],
                         ['A', 'C', None],
                         ['C', 'B', None]])

    print model.predict([['A', 'B', None],
                         ['A', None, 'C'],
                         [None, 'B', 'A']])

 

Relevant Link:  

https://yq.aliyun.com/articles/280788

  

7. 基於貝葉斯網進行圖像降噪

假設我們有以下圖像：

 

現在假設這張圖像受到了隨機噪聲的污染，變成了有噪聲的圖像：

現在我們的目標是恢復原始圖像。讓我們看看如何使用概率圖模型來實現。

0x1：問題分析

我們將有噪聲圖像中的每個像素都定義為一個觀察到的隨機變量，並將基准圖像中的每個像素都定義為一個未被觀察到的變量，由此，觀察到的變量和未被觀察到的變量都為 MxN 個。

我們將觀察到的變量表示為 X_ij，未被觀察到的變量定義為 Y_ij。每個變量都可取值 +1 或 -1（分別對應於黑色像素和白色像素）。

給定觀察到的變量，我們希望找到未觀察到的變量的最有可能的值。這對應於 MAP 推理。

0x2：圖結構構建

現在讓我們使用一些領域知識來構建圖結構：

• 在有噪聲圖像中的 (i,j) 位置觀察到的變量取決於在基准圖像中的 (i,j) 位置未觀察到的變量。原因是大多數時候它們是相等的，即 X_ij 和 Y_ij 的初始先驗是相等的。
• 對於基准圖像，
  • 在每個單一顏色的區域內，鄰近的像素通常有一樣的值，因此，如果 Y_ij 和 Y_kl 是鄰近像素，那么我們將它們連接起來。
  • 在顏色變化的邊界沒有這個性質。

由此，我們得到圖結構：

• 白色節點表示未被觀察到的變量 Y_ij
• 灰色節點表示觀察到的變量 X_ij
• 每個 X_ij 都連接到對應的 Y_ij
• 每個 Y_ij 都連接到它的相鄰節點

這個圖結構的先驗知識非常重要，這是我們后面進行概率推理的一個重要條件。

0x3：圖推理過程

我們的 MAP 推理問題可以用數學的方式寫出，如下： 

現在，我們需要根據我們的圖結構來定義我們的聯合分布 P(X,Y)。讓我們假設 P(X,Y) 由兩類因子組成，ϕ(X_ij, Y_ij) 和 ϕ(Y_ij,Y_kl)，對應於圖中的兩類邊。

接下來，我們按如下方式定義這些因子：

• ϕ(X_ij, Y_ij) = exp(w_e X_ij Y_ij)，其中 w_e 是一個大於 0 的參數。當 X_ij 和 Y_ij 相等時，這個因子取較大的值，當 X_ij 和 Y_ij 不同時就取較小的值。
• ϕ(Y_ij, Y_kl) = exp(w_s Y_ij Y_kl)，其中 w_s 也是一個大於 0 的參數。當 Y_ij 和 Y_kl 取值相等時，這個因子較大。

因此，我們的聯合分布可由下式給出：

將其插入到我們的 MAP 推理式中，可得：

有了這些參數后，我們需要求解上述 MAP 推理問題。我們可以使用迭代條件模式（ICM/iterated conditional mode）。其基本思想是：在每個步驟選擇一個節點 Y_ij，然后檢查 Y_ij=-1 和 Y_ij=1 時 MAP 推理表達式的值，然后選擇值更高的那個。重復這個過程一定的迭代次數或直到收斂，通常就能得到相當好的結果。

下面代碼中展示了這個過程，需要特別體會的是，這是一個動態的概率調整過程，每次像素的猜測的同時都相當於向網絡中輸入了新的信息，網絡中的條件概率也因此發生了改變。

# -*- coding: utf-8 -*-

import numpy as np
from PIL import Image
import sys


# look up the value of Y for the given indices
# if the indices are out of bounds, return 0
def compute_log_prob_helper(Y, i, j):
    try:
        return Y[i][j]
    except IndexError:
        return 0


def compute_log_prob(X, Y, i, j, w_e, w_s, y_val):
    result = w_e * X[i][j] * y_val

    # 在每個單一顏色的區域內，鄰近的像素通常有一樣的值，因此，如果 Y_ij 和 Y_kl 是鄰近像素，那么我們將它們連接起來
    result += w_s * y_val * compute_log_prob_helper(Y, i-1, j)
    result += w_s * y_val * compute_log_prob_helper(Y, i+1, j)
    result += w_s * y_val * compute_log_prob_helper(Y, i, j-1)
    result += w_s * y_val * compute_log_prob_helper(Y, i, j+1)
    return result


def denoise_image(X, w_e, w_s):
    m, n = np.shape(X)
    # initialize Y same as X
    # 在有噪聲圖像中的 (i,j) 位置觀察到的變量取決於在基准圖像中的 (i,j) 位置未觀察到的變量。原因是大多數時候它們是相等的，即 X_ij 和 Y_ij 的初始先驗是相等的
    Y = np.copy(X)
    # optimization
    max_iter = 10 * m * n
    for iter in range(max_iter):
        # randomly pick a location
        i = np.random.randint(m)
        j = np.random.randint(n)
        # compute the log probabilities of both values of Y_ij
        log_p_neg = compute_log_prob(X, Y, i, j, w_e, w_s, -1)
        log_p_pos = compute_log_prob(X, Y, i, j, w_e, w_s, 1)
        # assign Y_ij to the value with higher log probability
        if log_p_neg > log_p_pos:
            Y[i][j] = -1
        else:
            Y[i][j] = 1
        if iter % 100000 == 0:
          print ('Completed', iter, 'iterations out of', max_iter)
    return Y


# preprocessing step
def read_image_and_binarize(image_file):
    im = Image.open(image_file).convert("L")
    A = np.asarray(im).astype(int)
    A.flags.writeable = True

    A[A<128] = -1
    A[A>=128] = 1
    return A


def add_noise(orig):
    A = np.copy(orig)
    for i in range(np.size(A, 0)):
        for j in range(np.size(A, 1)):
            r = np.random.rand()
            if r < 0.1:
                A[i][j] = -A[i][j]
    return A


def convert_from_matrix_and_save(M, filename, display=False):
    M[M==-1] = 0
    M[M==1] = 255
    im = Image.fromarray(np.uint8(M))
    if display:
        im.show()
    im.save(filename)


def get_mismatched_percentage(orig_image, denoised_image):
    diff = abs(orig_image - denoised_image) / 2
    return (100.0 * np.sum(diff)) / np.size(orig_image)


def main():
    # read input and arguments
    orig_image = read_image_and_binarize("input.png")

    w_e = 8
    w_s = 10

    # add noise
    noisy_image = add_noise(orig_image)

    # use ICM for denoising
    denoised_image = denoise_image(noisy_image, w_e, w_s)

    # print the percentage of mismatched pixels
    print ('Percentage of mismatched pixels: ', get_mismatched_percentage(orig_image, denoised_image))

    convert_from_matrix_and_save(orig_image, 'orig_image.png', display=False)
    convert_from_matrix_and_save(noisy_image, 'noisy_image.png', display=False)
    convert_from_matrix_and_save(denoised_image, 'denoised_image.png', display=False)


if __name__ == '__main__':
    main()

 

降噪后的圖像

可以看到，算法獲得了70%的降噪效果。

Relevant Link:   

https://www.coursera.org/learn/probabilistic-graphical-models?authMode=signup
https://www.cnblogs.com/tornadomeet/p/3480694.html 
https://yq.aliyun.com/articles/280788  
https://github.com/prasoongoyal/image-denoising-mrf/blob/master/denoise.py