演算法實驗分析

Ⅰ EM演算法深度解析

最近在做文本挖掘的時候遇到了EM演算法，雖然讀書的時候簡單地接觸過，但當時並沒有深入地去了解，導致現在只記得演算法的名字。既然EM演算法被列為數據挖掘的十大演算法之一，正好借這個機會，重新學習一下這個經典的演算法。學習的過程中，我發現網上的資料大多講解地不夠細致，很多地方解釋得並不明了。因此我決定拋開別人的想法，僅從數學推導本身出發，盡力理解每一個公式的含義，並將其對應到實際的實驗過程當中。這篇博客記錄了我對與EM演算法的思考與理解，也是我人生中的第一篇博客，希望能夠對於想要學習EM演算法的同學有所幫助。

前面談到我在做文本挖掘的時候遇到了EM演算法，EM演算法用於估計模型中的參數。提到參數估計，最常見的方法莫過於極大似然估計——在所有的候選參數中，我們選擇的參數應該讓樣本出現的概率最大。相信看到這篇筆記的同學一定對極大似然估計非常熟悉，而EM演算法可以看作是極大似然估計的一個擴充，這里就讓我們用極大似然估計來解決一個簡單的例子，來開始正式的討論。

有A，B，C三枚硬幣，我們想要估計A，B，C三枚硬幣拋出正面的概率 , , 。我們按如下流程進行實驗100次：

記錄100次實驗的結果如下：

我們將上面的實驗結果表述如下：
表示第i次實驗中，硬幣A的結果，1代表正面，0代表反面； 表示第i次實驗中，硬幣B或硬幣C拋出正面的個數，則參數 的極大似然估計分別為：

即硬幣A，B，C各自拋出正面的次數占總次數的比例，其中 為指示函數。

實驗流程與1相同，但是我們不慎遺失了硬幣A的記錄結果，導致我們只知道隨後十次拋出了多少次正面，多少次反面，卻不知道實驗結果來自於硬幣B還是硬幣C。在這種情況下，我們是否還能估計出 , , 的值呢？

這時候利用極大似然估計似乎行不通了， 因為這種情況下，我們不但缺失了硬幣A產生的觀測值，同時也不知道哪些觀測值屬於硬幣B，哪些觀測值屬於硬幣C。

有些同學可能會提出，雖然我們無法得到三個硬幣各自產生的樣本，但是我們依然可以得到每個觀測值出現的概率。比如在第一次實驗中， 我們拋出了5次正面5次反面，我們可以做如下思考：

  假設這5次正面由硬幣B得到，那麼概率應該為 ，而這次觀測值來自於硬幣B，也就是硬幣A拋出正面的概率為

  假設這5次正面由硬幣C得到，那麼概率應該為 ，而這次觀測值來自於硬幣C，也就是硬幣A拋出反面的概率為

  綜合起來，利用條件概率公式，這個觀測值出現的概率就是

因此我們可以將樣本整體的概率和似然函數利用 , , 表示出來，通過對似然函數求導，令其關於 的偏導數等於0，我們可以求出三個參數的值。

這個思路聽上去十分合理，我們可以順著這個思路進行數學推導，看看可以得到什麼樣的結果。首先我們計算樣本的概率:

對應的似然函數為

其中 關於 的條件分布為

的分布為

因此我們可以得到

至此，我們成功地得到了似然函數。然而觀察可以發現，這個函數是由100項對數函數相加組成，每個對數函數內部包含一個求和，想通過求導並解出導數的零點幾乎是不可能的。當然我們可以通過梯度下降來極小化這個函數，藉助深度學習庫的自動微分系統在實現上也非常容易。但是這種做法過於簡單粗暴，有沒有辦法來優雅地解決這個問題呢？在繼續討論之前，我們先將這類問題進行一般化表述：

我們觀測到隨機變數 產生的m個相互獨立的樣本 , 的分布由聯合分布 決定， 是缺失數據或無法在實驗中被直接觀測到，稱為 隱變數 ，我們想要從樣本中估計出模型參數 的值。在接下來的討論中，我們假定 的取值是離散的，於是可以得到似然函數如下:

接下來，我們就探討一下，如何利用EM演算法解決這個問題。

這一部分的數學推導，主要參考了吳恩達CS229n的筆記，並且根據個人的思考和理解，盡力對公式的每一步進行詳細的解釋。我們先簡單地介紹一下琴生不等式。

琴生不等式有多種形式，下面給出其離散形式的表述和概率論中的表述:
1.若 為嚴格凹函數， 為定義域內的n個點， 是n個正實數，且滿足 , 則下述不等式成立:

當且僅當 時，不等式取等號。

2.若 為嚴格凹函數， 為實值隨機變數，且期望存在，則下述不等式成立:

當且僅當 ，即 為常數時，不等式取等號。

註： 這里將函數上方為凹集的函數稱為凹函數， 例如 函數就是凹函數。
相信大家對琴生不等式都十分熟悉，因此這里就不做過多的說明。接下來，我們將琴生不等式應用到我們的問題中。

回到我們之前的問題上， 我們想要極大化下面這個函數:

但是我們無法對這個函數直接求導，因此我們藉助琴生不等式，對這個函數進行變換。為了讓過程看上去簡潔，下面只對求和中的第 項進行計算。

令 滿足 ，且 ，則根據琴生不等式，可以得到：

當且僅當 為常數時，上述不等式取等號。也就是說，對於任意 ， 是一個與 無關的量。設對於任意 ，我們可以得到：

因此當 時，不等式 取等號，容易驗證此時 , 與 無關。將 綜合一下，我們可以得到以下結論:

到這里為止，我們已經擁有了推導出EM演算法的全部數學基礎，基於 我們可以構建出E步和M步。上面的數學推導雖然看上去略為復雜，但實際上只用到了三個知識點：
  1.琴生不等式:

  2.條件概率:

  3.聯合分布求和等於邊緣分布:

對上面的數學推導有疑問的同學，可以結合上面這三點，再將整個推導過程耐心地看一遍。

大部分關於EM演算法的資料，只是在數學形式上引入了 函數，即 ，以滿足琴生不等式的使用條件，卻沒有過多地解釋 函數本身。這導致了很多人完全看懂了演算法的推導，卻還是不理解這些數學公式究竟在做什麼，甚至不明白EM演算法為什麼叫做EM演算法。所以在給出E步和M步之前，我想先談一談 函數。

我們回顧一下 函數所滿足的條件（暫時不考慮琴生不等式取等號的限制），

在 所有可能的取值處有定義。可以看出， 是 的樣本空間上任意的一個概率分布。因此，我們可以對不等式 進行改寫。首先我們可以將含有 的求和寫成期望的形式:

這里 指的是在概率分布 下，求隨機變數 和 的期望。有同學會問，為什麼我們平時求期望的時候只要寫 ，並沒有指明是在哪個概率分布下的期望。這是因為一般情況下，我們都清楚地知道隨機變數 所服從的分布 ，並且默認在分布 下求期望。

舉個例子，我手上有一個硬幣，拋了10次，問拋出正面次數的期望。這種情況下，大部分人會默認硬幣是均勻的，也就是說拋出正面的次數 服從二項分布 ，期望 。這時有人提出了質疑，他說我認為你這個硬幣有問題，拋出正面的概率只有0.3，那麼在他眼裡， 期望 。

回到正題，我們利用等式 改寫不等式 ，可以得到:

這正是琴生不等式在概率論中的形式。我們可以將不等式倒過來理解：
  首先，假定隨機變數 服從概率分布 ， 是 的樣本空間上的任意一個概率分布。這里 可以是一組定值，也可以是關於參數 的函數。

  顯然，當我們取不同的 時，隨機變數 的期望也會隨之改變。需要注意的是，由於 與 相關，所以這里的期望不是一個數值，而是關於 的函數。

  當我們令 為 的後驗分布 時，上面的期望最大。這里有兩點需要注意，1. 後驗分布 也是一個關於參數 的函數。2. 由於期望是關於 的函數，所以這里的最大指的並非是最大值，而是最大的函數。

  若對於每一個 ，我們都令 為 的後驗分布 ，則上述期望之和等於我們要極大化的似然函數，即

通過上述分析，我們為尋找似然函數的極大值點 提供了一個思路。我們不去極大化似然函數本身，而是去極大化 。至於如何將這個思路實際應用，就要利用到EM演算法中的E-step和M-step。

這一節中，我們先給出E-step和M-step的數學形式，隨後在結合拋硬幣的例子來解釋這兩步究竟在做什麼。下面進入演算法的流程，首先我們任意初始化 ，按下述過程進行迭代直至收斂：

在第 次迭代中，
(E-step)對於每個 ，令
(M-step)更新 的估計值，令

EM演算法從任意一點 出發，依次利用E-step優化 ，M-step優化 ，重復上述過程從而逐漸逼近極大值點。而這個過程究竟是怎樣的呢，就讓我們一步步地揭開EM演算法的面紗。

假設我們現在隨機初始化了 ，進入第一輪迭代：
(E-step)

由於我們已經假定模型參數為 ，所以此時 不再是與 有關的函數，而是由一組常數構成的概率分布。結合拋硬幣的例子來看，這一步是在我們已知模型參數 的基礎上(雖然這是我們瞎猜的)，去推測每一次的觀測值是由哪個硬幣產生的，或者說我們對每一次觀測值做一個軟分類。比如我們根據初始化的參數，計算出 ， 。可以解釋為第 個觀測值有20%的概率來自於硬幣B，80%的概率來自於硬幣C；或者說硬幣A拋出了0.2個正面，0.8個反面。

(M-step)

考慮到 是一組常數，我們可以舍棄常數項，進一步簡化上面這個要極大化的函數

由於 不再與 相關，因此上面的函數變成了對數函數求和的形式，這個函數通常來說是容易求導的，令導數等於0，我們可以求出新的參數 。我們仍舊以拋硬幣為例進行解釋，

令 , 可以得到，

這三個參數的解釋是顯而易見的。我們在E-step中對每個觀測值進行了軟分類， 可以看成是硬幣A拋出正面的次數，所以 是 的極大似然估計； 是我們拋硬幣B的次數， 是硬幣B拋出正面的次數，所以 是 的極大似然估計；對於 我們有相同的解釋。

我們將這個結果與拋硬幣1中極大似然估計的結果相比較可以發現，之前結果中的指示函數 變成了這里的 ，在指示函數下，某個觀測值要麼來自於硬幣B，要麼來自於硬幣C，因此也稱為硬分類。而在 函數下，某個觀測值可以一部分來自於硬幣B，一部分來自於硬幣C，因此也稱作軟分類。

將上述兩步綜合起來，EM演算法可以總結如下：我們首先初始化模型的參數，我們基於這個參數對每一個隱變數進行分類，此時相當於我們觀測到了隱變數。有了隱變數的觀測值之後，原來含有隱變數的模型變成了不含隱變數的模型，因此我們可以直接使用極大似然估計來更新模型的參數，再基於新的參數開始新一輪的迭代，直到參數收斂。接來下我們就討論為什麼參數一定會收斂。

前面寫了太多的公式，但是這一部分我不打算給出收斂性的數學推導。其實數學上證明EM演算法的收斂性很容易，只需要證明每一輪迭代之後，參數的似然函數遞增，即

Ⅱ 07 - 密度聚類演算法（DBSCAN）實驗案例

DBSCAN演算法實驗案例分析

DBSCAN，作為一種基於密度的聚類演算法，以其在發現樣本密集區域和處理任意形狀數據集的能力而著稱。本實驗通過三個實際案例來展示其效能：環形數據、新月形數據以及輪廓系數的評估。

環形數據聚類

首先，通過scikit-learn生成環形數據，目標是分成3類。K-Means、MeanShift和Birch演算法未能達到預期，而未調參的DBSCAN將所有數據歸為一類。通過調整eps和min_samples參數，DBSCAN成功實現了分類。

新月形數據聚類

新月形數據集需要分為上下兩個類別。在嘗試多種演算法後，DBSCAN憑借其靈活性，通過參數調整，成功完成了任務，體現了對任意形狀數據的適應性。

輪廓系數評估

使用輪廓系數對聚類效果進行量化評估，這個指標關注樣本的內聚度和分離度，值越高表示聚類效果越好。通過計算，可以客觀地對比不同演算法的聚類效果。

Ⅲ 貪心演算法的例題分析

例題1、
[0-1背包問題]有一個背包，背包容量是M=150。有7個物品，物品不可以分割成任意大小。
要求盡可能讓裝入背包中的物品總價值最大，但不能超過總容量。
物品 A B C D E F G
重量 35kg 30kg 6kg 50kg 40kg 10kg 25kg
價值 10$ 40$ 30$ 50$ 35$ 40$ 30$
分析：
目標函數：∑pi最大
約束條件是裝入的物品總重量不超過背包容量：∑wi<=M(M=150）
⑴根據貪心的策略，每次挑選價值最大的物品裝入背包，得到的結果是否最優？
⑵每次挑選所佔重量最小的物品裝入是否能得到最優解？
⑶每次選取單位重量價值最大的物品，成為解本題的策略。
值得注意的是，貪心演算法並不是完全不可以使用，貪心策略一旦經過證明成立後，它就是一種高效的演算法。
貪心演算法還是很常見的演算法之一，這是由於它簡單易行，構造貪心策略不是很困難。
可惜的是，它需要證明後才能真正運用到題目的演算法中。
一般來說，貪心演算法的證明圍繞著：整個問題的最優解一定由在貪心策略中存在的子問題的最優解得來的。
對於例題中的3種貪心策略，都是無法成立（無法被證明）的，解釋如下：
⑴貪心策略：選取價值最大者。
反例：
W=30
物品：A B C
重量：28 12 12
價值：30 20 20
根據策略，首先選取物品A，接下來就無法再選取了，可是，選取B、C則更好。
⑵貪心策略：選取重量最小。它的反例與第一種策略的反例差不多。
⑶貪心策略：選取單位重量價值最大的物品。
反例：
W=30
物品：A B C
重量：28 20 10
價值：28 20 10
根據策略，三種物品單位重量價值一樣，程序無法依據現有策略作出判斷，如果選擇A，則答案錯誤。
【注意：如果物品可以分割為任意大小，那麼策略3可得最優解】
對於選取單位重量價值最大的物品這個策略，可以再加一條優化的規則：對於單位重量價值一樣的，則優先選擇重量小的！這樣，上面的反例就解決了。
但是，如果題目是如下所示，這個策略就也不行了。
W=40
物品：A B C
重量：25 20 15
價值：25 20 15
附：本題是個DP問題，用貪心法並不一定可以求得最優解，以後了解了動態規劃演算法後本題就有了新的解法。
例題2、
馬踏棋盤的貪心演算法
123041-23 XX
【問題描述】
馬的遍歷問題。在8×8方格的棋盤上，從任意指定方格出發，為馬尋找一條走遍棋盤每一格並且只經過一次的一條路徑。
【初步設計】
首先這是一個搜索問題，運用深度優先搜索進行求解。演算法如下：
⒈ 輸入初始位置坐標x,y；
⒉ 步驟 c:
如果c> 64輸出一個解，返回上一步驟c--
（x,y) ← c
計算(x,y）的八個方位的子結點，選出那些可行的子結點
循環遍歷所有可行子結點，步驟c++重復2
顯然⑵是一個遞歸調用的過程，大致如下：
C++程序： #defineN8voiddfs(intx,inty,intcount){inti,tx,ty;if(count>N*N){output_solution();//輸出一個解return;}for(i=0;i<8;i++){tx=hn[i].x;//hn[]保存八個方位子結點ty=hn[i].y;s[tx][ty]=count;dfs(tx,ty,count+1);//遞歸調用s[tx][ty]=0;}}Pascal程序： ProgramYS;ConstFXx:array[1..8]of-2..2=（1,2,2,1,-1,-2,-2,-1）；FXy:array[1..8]of-2..2=（2,1,-1,-2,-2,-1,1,2）；VarRoad:array[1..10,1..10]ofinteger;x,y,x1,y1,total:integer;ProcereFind(x,y:integer);varNx,Ny,i:integer;BeginFori:=1to8dobegin{8個方向}If(x+FXx[i]in[1..8])and(y+FXy[i]in[1..8])Then{確定新坐標是否越界}IfRoad[x+Fxx[i],y+Fxy[i]]=0Thenbegin{判斷是否走過}Nx:=x+FXx[i];Ny:=y+FXy[i];Road[Nx,Ny]:=1;{建立新坐標}If(Nx=x1）and(Ny=y1）Theninc(total)elseFind(Nx,Ny);{遞歸}Road[Nx,Ny]:=0{回朔}endendEnd;BEGIN{Main}Total:=0;FillChar(Road,sizeof(road),0);Readln(x,y);{讀入開始坐標}Readln(x1,y1）；{讀入結束坐標}If(x>10）or(y>10）or(x1>10）or(y1>10）Thenwriteln('Error'){判斷是否越界}ElseFind(x,y);Writeln('Total:',total){打出總數}END.這樣做是完全可行的，它輸入的是全部解，但是馬遍歷當8×8時解是非常之多的，用天文數字形容也不為過，這樣一來求解的過程就非常慢，並且出一個解也非常慢。
怎麼才能快速地得到部分解呢？
【貪心演算法】
其實馬踏棋盤的問題很早就有人提出，且早在1823年，J.C.Warnsdorff就提出了一個有名的演算法。在每個結點對其子結點進行選取時，優先選擇『出口』最小的進行搜索，『出口』的意思是在這些子結點中它們的可行子結點的個數，也就是『孫子』結點越少的越優先跳，為什麼要這樣選取，這是一種局部調整最優的做法，如果優先選擇出口多的子結點，那出口少的子結點就會越來越多，很可能出現『死』結點（顧名思義就是沒有出口又沒有跳過的結點），這樣對下面的搜索純粹是徒勞，這樣會浪費很多無用的時間，反過來如果每次都優先選擇出口少的結點跳，那出口少的結點就會越來越少，這樣跳成功的機會就更大一些。這種演算法稱為為貪心演算法，也叫貪婪演算法或啟發式演算法，它對整個求解過程的局部做最優調整，它只適用於求較優解或者部分解，而不能求最優解。這樣的調整方法叫貪心策略，至於什麼問題需要什麼樣的貪心策略是不確定的，具體問題具體分析。實驗可以證明馬遍歷問題在運用到了上面的貪心策略之後求解速率有非常明顯的提高，如果只要求出一個解甚至不用回溯就可以完成，因為在這個演算法提出的時候世界上還沒有計算機，這種方法完全可以用手工求出解來，其效率可想而知。